ASER Pakistan 2016

In this piece of paper, a set of data obtained from Annual Status of Education Report (ASER) is explored. The raw data was downloaded from the link here. https://palnetwork.org/aser-centre/

Preparation

Packages Used

library(tidyverse)
library(ggplot2)
library(ggthemes)
library(ggrepel)
library(gghighlight)
library(stringr)
library(dplyr)
library(sf)
library(scatterplot3d)
library(car)
library(ResourceSelection) #to excute Hosmer-Lemeshow test
## Warning: package 'ResourceSelection' was built under R version 4.0.3
# library(broom)
require(ggiraph)
## Warning: package 'ggiraph' was built under R version 4.0.3
require(ggiraphExtra)
# require(plyr)

Data Installation

provdist <- read.csv("aser/ASER2016ProvDist.csv")
school <- read.csv("aser/ASER2016GSchool.csv")
child <- read.csv("aser/ASER2016Child.csv")
pschool <- read.csv("aser/ASER2016PvtSchool.csv")
gschool <- read.csv("aser/ASER2016GSchool.csv")
parent <- read.csv("aser/ASER2016Parent.csv")
house <- read.csv("aser/ASER2016HouseholdIndicators.csv")
RegionName <- c("2" = "Panjab", 
                "3" = "Sindh", 
                "4" = "Balochistan", 
                "5" = "Khyber Pakhtunkhwa", 
                "6" = "Gilgit-Baltistan", 
                "7" = "Azad Jammu and Kashmir", 
                "8" = "Islamabad - ICT", 
                "9" = "Federally Administrated Tribal Areas")
Gender <- c("0" = "Male",
            "-1" = "Female")

Exploration

Checking Samplesizes

length(unique(child$CID))
## [1] 255196

The whole samplesize (the numebr of children) of this dataset is 255196.

child %>% 
  filter(DID == 266) %>% 
  summarize(N_hunza = length(unique(CID)))

The samplesize of Hunza alone is 1641.

Exploration in Hunza

Gender Proportion

child %>% 
  filter(DID == 266) %>% 
  summarize(gender_proportion = mean(C002))

-1: female, 0: male
gender_proportion = -0.5173675 means there are a little more girls in the dataset.

Age of Children (C001)

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C001)) +
  geom_histogram()

Age is well sparsed

Eduation Status

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C003)) +
  geom_histogram(bins = 3)

1 = never enrolled; 2 = drop-out; 3 = currently enrolled

Education Status by Gender

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C003)) +
  geom_histogram(bins = 3, binwidth = 1) +
  facet_grid(~C002, labeller = labeller(C002 = Gender))

Both genders look pretty good interms of the absolute number of currently-enrolled-children

The Enrollment Rate by Gender

child %>% 
  filter(DID == 266) %>% 
  group_by(C002) %>% 
  mutate(enrollment_rate = mean(C003 == 3)) %>% 
  ggplot(aes(C002, enrollment_rate)) +
  geom_col() +
  scale_y_continuous() +
  geom_label(aes(label = enrollment_rate)) +
  scale_x_continuous(breaks = c(-1, 0), labels = c("Female", "Male"))

As a rate, both are doing pretty good

Currently-Enrolled: Institution Type (C006) in Hunza

1 = Government; 2 = Private; 3 = Madrasah(Conventional religious education) School; 4 = Other(Non formal education facility)

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C006)) +
  geom_histogram()

Most children go to public schools or private schools

Within Gilgit-Baltistan

child %>% 
  filter(RID == 6) %>% 
  group_by(DID) %>% 
  mutate(Current_Enrollment_Rate = mean(C003 == 3)) %>% 
  ggplot(aes(DID, Current_Enrollment_Rate)) +
  geom_count() +
  scale_x_continuous(breaks = 260:266, labels = c("Gilgit", "Diamer", "Skardu", "Ghanshe", "Astore", "Ghizer", "Hunza-Nagar"))

Within Gilgit-Baltistan, Hunza is outperforming.

Basic Learning Levels: Reading in Local/National Language (C010)

1 = Begginer/Nothing; 2 = Letters; 3 = Words; 4 = Sentences; 5 = Story

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C010)) +
  geom_histogram()

child %>% 
  filter(DID == 266) %>% 
  summarize(na = sum(is.na(C010)))

Basic Learning Levels by Age

child %>% 
  filter(DID == 266, C013 != c(3,4)) %>% 
  ggplot(aes(C010)) +
  geom_histogram() +
  facet_grid(~C001)

# children at he age of 3 and 4 are removed for they have not data

Basic Learning Levels by Gender

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C010)) +
  geom_histogram() +
  facet_grid(~C002, labeller = labeller(C002 = Gender))

English Reading Levels (C013)

child %>% 
  filter(DID == 266, C013 != c(3,4)) %>% 
  ggplot(aes(C013)) +
  geom_histogram() +
  facet_grid(~C001)

# children at he age of 3 and 4 are removed for they have not data

English Reading Levels by Gender

child %>% 
  filter(DID == 266) %>% 
  ggplot(aes(C013)) +
  geom_histogram() +
  facet_grid(~C002, labeller = labeller(C002 = Gender))

Comparison between Other Region

Current Enrollment Rate

child %>% 
  group_by(DID) %>% 
  mutate(avg = round(mean(C003 == 3), digits = 2)) %>% 
  ungroup() %>% 
  ggplot(aes(avg)) +
  geom_histogram() +
  facet_grid(~RID, labeller = labeller(RID = RegionName)) +
  labs(title = "Current Enrollment Rate by Region")

Female Average Learning Levels (Age or other variables are NOT adjusted)

child %>% 
  filter(C002 == -1) %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(DID, avg_learning, color = RID)) +
  geom_point() +
  geom_text(aes(label = DID), nudge_x = 5, check_overlap = TRUE)

child %>% 
  filter(C002 == -1) %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(DID, avg_learning, color = RID)) +
  geom_point() +
  geom_text(aes(label = DID), nudge_x = 5, check_overlap = TRUE) +
  gghighlight(RID == 6)

It is interesting to note that Gilgit-Baltistan(RID==6) has a huge diversity in average learning levels of girls and Hunza(DID==266) is in the top group of all region.

Spatial Analysis

Districts Map

ica_df <- ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
  as.data.frame()
## Warning: Problem with `mutate()` input `centroid`.
## x st_centroid does not give correct centroids for longitude/latitude data
## i Input `centroid` is `st_centroid(geometry)`.
## Warning in st_centroid.sfc(geometry): st_centroid does not give correct
## centroids for longitude/latitude data
ica_df <- ica_df %>% select(Province, Districts, x, y)
ica_df <- ica_df %>% summarize(Province = tolower(Province), Districts = tolower(Districts), x = x, y = y)

Child and ProvDist data combination

child_dname <- child %>% left_join(provdist[-1])
## Joining, by = "DID"
child_dname <- child_dname %>% mutate(dname = tolower(DNAME))
ica_df_3 <- ica_df %>% filter(Province == "sindh")

ica_df_3$Districts <- ica_df_3$Districts %>%  
  str_replace("ghotki", "gotki") %>%
  str_replace("mirpur khas", "mirpurkhas") %>% 
  str_replace("malir karachi", "karachi-malir-rural") %>% 
  str_replace("naushahro feroze", "nowshero feroze") %>% 
  str_replace("kambar shahdad kot", "qambar shahdadkot") %>% 
  str_replace("sujawal", "sajawal") %>% 
  str_replace("shaheed benazir abad", "shaheed benazirabad") %>% 
  str_replace("tando allahyar", "tando allah yar") %>% 
  as.vector()

child_dname_3 <- child_dname %>% filter(RNAME == "Sindh") %>% left_join(ica_df_3, by = c("dname" = "Districts"))

child_dname_3 %>% group_by(dname) %>% summarize(n = sum(x))
## `summarise()` ungrouping output (override with `.groups` argument)
ica_df_3

District IDs

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica, aes(x, y, label = DID, color = factor(RID))) +
  geom_text(data = child_ica, aes(x, y, label = DID), size = 5, nudge_y = 0.2, check_overlap = TRUE)
## Warning: Problem with `mutate()` input `centroid`.
## x st_centroid does not give correct centroids for longitude/latitude data
## i Input `centroid` is `st_centroid(geometry)`.
## Warning in st_centroid.sfc(geometry): st_centroid does not give correct
## centroids for longitude/latitude data
## Warning: Ignoring unknown aesthetics: label
## Warning: Removed 8541 rows containing missing values (geom_point).
## Warning: Removed 8541 rows containing missing values (geom_text).

Gilgit-Baltistan Region

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica, aes(x, y, label = DNAME, shape = factor(RID), color = factor(RID))) +
  geom_text(data = child_ica, aes(x, y, label = DNAME), size = 4.5, nudge_y = 0.2, check_overlap = TRUE) +
  scale_shape_manual(values = 0:8) +
  coord_sf(xlim = c(71, 77.9), ylim = c(34, 37.1))
## Warning: Problem with `mutate()` input `centroid`.
## x st_centroid does not give correct centroids for longitude/latitude data
## i Input `centroid` is `st_centroid(geometry)`.
## Warning in st_centroid.sfc(geometry): st_centroid does not give correct
## centroids for longitude/latitude data
## Warning: Ignoring unknown aesthetics: label
## Warning: Removed 8541 rows containing missing values (geom_point).
## Warning: Removed 8541 rows containing missing values (geom_text).

Latitudes and longitudes of Gilgit-Baltistan region

child_ica %>% 
  filter(RID == 6) %>% 
  summarise(xmin = min(x), xmax = max(x), ymin = min(y), ymax = max(y))

Average enrollment rate

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica %>% group_by(DID) %>% 
               mutate(avg_enroll = mean(C003 == 3)), 
             aes(x, y, label = avg_enroll, color = avg_enroll))

  # geom_text(data = child_ica%>% group_by(DID) %>% 
  #             mutate(avg_enroll = mean(C003 == 3)), 
  #           aes(x, y, avg_enroll), check_overlap = TRUE, nudge_y = 0.5)

Average enrollment rate of female children

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica %>% filter(C002 == -1) %>% group_by(DID) %>% 
               mutate(avg_enroll = mean(C003 == 3)), 
             aes(x, y, label = avg_enroll, color = avg_enroll))
## Warning: Problem with `mutate()` input `centroid`.
## x st_centroid does not give correct centroids for longitude/latitude data
## i Input `centroid` is `st_centroid(geometry)`.
## Warning in st_centroid.sfc(geometry): st_centroid does not give correct
## centroids for longitude/latitude data
## Warning: Ignoring unknown aesthetics: label
## Warning: Removed 3566 rows containing missing values (geom_point).

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica, aes(x, y, label = C010, color = C010))

  # geom_text(data = child_ica, aes(x, y, label = C010), check_overlap = TRUE, nudge_y = 1)

Traial: children, gender

ica %>% 
  mutate(centroid = st_centroid(geometry),
    x = st_coordinates(centroid)[,1],
    y = st_coordinates(centroid)[,2]) %>% 
    ggplot() +
  geom_sf() +
  geom_point(data = child_ica %>% 
               group_by(DID) %>%
               mutate(gender_ratio = mean(C002)), 
             aes(x, y, label = gender_ratio, color = gender_ratio))

0: male, -1: female

Thus, -0.5 indicates the complete gender parity

Children

Average Local/National Language Level Across Ages by Rgion

child_ica %>% 
  filter(C002 == c(0,-1)) %>% 
  group_by(DID, C001, C002) %>% 
  mutate(avg_local = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_local, color = factor(C002), group = DID)) +
  geom_line(show.legend = FALSE) +
  labs(x = "Age") +
  facet_grid(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  facet_wrap(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender))

Average Enrollment Rate

Across Ages by Rgion
child_ica %>% 
  filter(C002 == c(0,-1)) %>% 
  group_by(DID, C001, C002) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, color = factor(DID), group = DID)) +
  geom_line(show.legend = FALSE) +
  labs(x = "Age") +
  facet_grid(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  facet_wrap(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender))

Hunza
child_ica %>% 
  filter(C002 != "NA", RID == 6, RID != "NA") %>% 
  group_by(DID, C001, C002) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, color = factor(DID), group = DID)) +
  geom_line() +
  labs(x = "Age") +
  facet_grid(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  facet_wrap(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  gghighlight(label_key = DNAME, calculate_per_facet = TRUE)

child_ica %>% 
  filter(C002 != "NA", RID == 6, RID != "NA") %>% 
  group_by(DID, C001, C002) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, color = factor(DID), group = DID)) +
  geom_line() +
  labs(x = "Age") +
  facet_grid(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  facet_wrap(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  gghighlight(DID == 266,label_key = DNAME, calculate_per_facet = TRUE)
## Warning: Tried to calculate with group_by(), but the calculation failed.
## Falling back to ungrouped filter operation...

Karachi
child_ica %>% 
  filter(C002 != "NA", RID == 3) %>% 
  group_by(DID, C001, C002) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, color = factor(DID), group = DID)) +
  geom_line() +
  labs(x = "Age") +
  facet_grid(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  facet_wrap(RID~C002, labeller = labeller(RID = RegionName, C002 = Gender)) +
  gghighlight(DID == c(315,316), label_key = DNAME, calculate_per_facet = TRUE)

children Education Stage
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", C010 != "NA", C002 != "NA") %>% 
  group_by(C005) %>% 
  summarize(n = n()) %>%
  mutate(C005 = fct_reorder(C005, n)) %>% 
  ggplot(aes(C005, n)) +
  geom_col() +
  coord_flip()
## `summarise()` ungrouping output (override with `.groups` argument)

Hunza. Children Education Stage
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", !is.na(C010), C002 != "NA", !is.na(C005), DID == 266) %>% 
  group_by(C005) %>% 
  summarize(n = n()) %>%
  mutate(C005 = fct_reorder(C005, n)) %>% 
  ggplot(aes(C005, n)) +
  geom_col() +
  coord_flip()
## `summarise()` ungrouping output (override with `.groups` argument)

Enrollment and Gender within Gilgit-Baltistan
child_ica %>% 
  filter(!is.na(C002), RID == 6) %>% 
  group_by(DID, C001, C002) %>% 
  mutate(enroll = mean(C003 == 3)) %>% 
  ggplot(aes(C001, enroll, group = C002, color = factor(C002))) +
  geom_line() +
  facet_wrap(.~DID)

Education Etages
child_ica %>% 
  group_by(DID) %>% 
  mutate(n = length(C003),
            dropout = sum(C003 == 2),
            dropout_ratio = mean(C003 == 2),
            never = sum(C003 == 1),
            never_ratio = mean(C003 == 1)) %>% 
  ggplot(aes(DID, never_ratio, color = factor(RID))) +
  geom_point()

Child & Parents

Mother Enrollment by District
child_ica %>%
  filter(RID == 6, PR004 != "NA") %>% 
  group_by(DID, PR004) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(DID, avg_enrollment, fill = factor(PR004))) +
  geom_col(aes(factor(DID), avg_enrollment, fill = factor(PR004)), position = "dodge")

Within Gilgit-Baltistan, Hunza is the only district where there are more women who have experience of education than those who have never enrolled in schools.

Father Enrollment by District
child_ica %>%
  filter(RID == 6, PR009 != "NA") %>% 
  group_by(DID, PR009) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(DID, avg_enrollment, group = PR009, fill = factor(PR009))) +
  geom_col(aes(factor(DID), avg_enrollment), position = "dodge")

Mother’s Education and Children Learning Level by Region. boxplot
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3)) %>% 
  ggplot(aes(PR004, avg_learning, group = PR004)) +
  geom_boxplot(aes(PR004, avg_learning)) +
  facet_grid(.~RID)

Father’s Education and Children Learning Level by Region. boxplot
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3)) %>% 
  ggplot(aes(PR009, avg_learning, group = PR009)) +
  geom_boxplot(aes(PR009, avg_learning)) +
  facet_grid(.~RID)

Mother’s Education and Children Learning Level by District
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3, na.rm = TRUE), mother_enrollment = mean(PR004, na.rm = TRUE)) %>% 
  ggplot(aes(mother_enrollment, avg_learning)) +
  geom_point()

child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3, na.rm = TRUE), mother_enrollment = mean(PR004, na.rm = TRUE)) %>% 
  ungroup() %>% 
  summarize(r = cor(mother_enrollment, avg_learning))
Father’s Education and Children Learning Level by District
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3, na.rm = TRUE), father_enrollment = mean(PR009, na.rm = TRUE)) %>% 
  ggplot(aes(father_enrollment, avg_learning)) +
  geom_point()

child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA") %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C003 == 3), father_enrollment = mean(PR009)) %>% 
  ungroup() %>% 
  summarize(r = cor(father_enrollment, avg_learning))
Hunza. Age, Learning Level, Mother’s Education,
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA", C003 != "NA", DID == 266, C001 != "NA") %>% 
  group_by(PR004, C001) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_learning, group = PR004, color = factor(PR004))) +
  geom_line(aes(C001, avg_learning))
## Warning: Removed 200 row(s) containing missing values (geom_path).

Hunza. Age, Learning Level, Father’s Education
child_ica %>% 
  filter(PR009 != "NA", PR004 != "NA", C003 != "NA", DID == 266, C001 != "NA") %>% 
  group_by(PR009, C001) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_learning, group = PR009, color = factor(PR009))) +
  geom_line(aes(C001, avg_learning))
## Warning: Removed 200 row(s) containing missing values (geom_path).

Gilgit-Baltistan. Age, Learning Level, Mother’s Education
child_ica %>% 
  filter(PR009 != "NA", 
         PR004 != "NA", 
         C003 != "NA", 
         C001 != "NA",
         RID == 6) %>% 
  group_by(DID, PR004, C001) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_learning, group = PR004, color = factor(PR004))) +
  geom_line(aes(C001, avg_learning)) +
  facet_grid(.~DID) +
  facet_wrap(.~DID)
## Warning: Removed 226 row(s) containing missing values (geom_path).

Gilgit-Baltistan. Age, Learning Level, Father’s Education
child_ica %>% 
  filter(PR009 != "NA", 
         PR004 != "NA", 
         C003 != "NA", 
         C001 != "NA",
         RID == 6) %>% 
  group_by(DID, PR009, C001) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_learning, group = PR009, color = factor(PR009))) +
  geom_line(aes(C001, avg_learning)) +
  facet_grid(.~DID) +
  facet_wrap(.~DID)
## Warning: Removed 226 row(s) containing missing values (geom_path).

Gilgit-Baltistan. Age, Enrollment, Mother’s Education
child_ica %>% 
  filter(PR009 != "NA", 
         PR004 != "NA", 
         C003 != "NA", 
         C001 != "NA",
         RID == 6) %>% 
  group_by(DID, PR004, C001) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, group = PR004, color = factor(PR004))) +
  geom_line(aes(C001, avg_enrollment)) +
  facet_grid(.~DID) +
  facet_wrap(.~DID)

Gilgit-Baltistan. age, enrollment, father’s education
child_ica %>%
  filter(PR009 != "NA", 
         PR004 != "NA", 
         C003 != "NA", 
         C001 != "NA",
         RID == 6) %>% 
  group_by(DID, PR009, C001) %>% 
  mutate(avg_enrollment = mean(C003 == 3, na.rm = TRUE)) %>% 
  ggplot(aes(C001, avg_enrollment, group = PR009, color = factor(PR009))) +
  geom_line(aes(C001, avg_enrollment)) +
  facet_grid(.~DID) +
  facet_wrap(.~DID)

mothers education level
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", C010 != "NA", C002 != "NA") %>% 
  ggplot(aes(PR005)) +
  geom_histogram(stat = "count") +
  coord_flip()
## Warning: Ignoring unknown parameters: binwidth, bins, pad

summary
child_ica %>% 
  ggplot(aes(PR006)) +
  geom_histogram(stat = "count") +
  coord_flip()
## Warning: Ignoring unknown parameters: binwidth, bins, pad

fathers education level
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", C010 != "NA", C002 != "NA") %>%
  ggplot(aes(PR010)) +
  geom_histogram(stat = "count") +
  coord_flip()
## Warning: Ignoring unknown parameters: binwidth, bins, pad

###### summary

child_ica %>% 
  ggplot(aes(PR011)) +
  geom_histogram(stat = "count") +
  coord_flip()
## Warning: Ignoring unknown parameters: binwidth, bins, pad

mothers education level by region
child_ica %>%
  group_by(DID, PR006) %>% 
  ggplot(aes(PR006)) +
  geom_histogram(aes(PR006), stat = "count") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, size = 5)) +
  facet_grid(.~RID)
## Warning: Ignoring unknown parameters: binwidth, bins, pad

Mother enrollment by district
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", C010 != "NA") %>% 
  group_by(DID) %>% 
  mutate(mother_edu_rate = sum(PR004 == -1)/length(PR004)) %>% 
  ggplot(aes(DID, mother_edu_rate, shape = factor(RID), color = factor(RID))) +
  geom_point() +
  scale_shape_manual(values = 0:8)

Gilgit-Baltistan. Average mother education experience, House type
child_ica %>% 
  filter(PR004 != "NA", PR009 != "NA", C010 != "NA") %>% 
  group_by(DID) %>% 
  mutate(mother_edu_rate = sum(PR004 == -1)/length(PR004)) %>% ungroup() %>% 
  ggplot(aes(factor(RID), mother_edu_rate, group = RID)) +
  geom_violin()

Hunza. father enroll ratio
child_ica %>% 
  filter(DID == 266, !is.na(PR009)) %>% 
  summarize(father_edu_ratio = sum(PR009 == -1)/length(PR009))
father enroll ratio
child_ica %>% 
  filter(!is.na(PR009)) %>% 
  summarize(father_edu_ratio = sum(PR009 == -1)/length(PR009))
father enroll ratio by district
child_ica %>% 
  filter(!is.na(PR009)) %>% 
  group_by(DID) %>% 
  mutate(father_edu_ratio = sum(PR009 == -1)/length(PR009)) %>% 
  ggplot(aes(factor(DID), father_edu_ratio)) +
  geom_point()

mother enroll ratio by district
child_ica %>% 
  filter(!is.na(PR009), !is.na(PR004)) %>% 
  group_by(DID) %>% 
  mutate(mother_edu_ratio = sum(PR004 == -1)/length(PR004)) %>% 
  ggplot(aes(factor(DID), mother_edu_ratio)) +
  geom_point()

Child, Household

Number of Children in a household
child_ica %>%
  filter(!is.na(PR004), !is.na(PR009)) %>% 
  filter(HHID == 96546) %>% 
  select(CID, PRID, C002, C001, PR001, VID, C003, DID, RID)
child_ica %>% 
  group_by(HHID) %>%
  mutate(n = n()) %>% 
  ggplot(aes(factor(n), fill = factor(RID))) +
  geom_bar(position = "dodge") +
  facet_grid(.~RID) +
  facet_wrap(.~RID)

house type ratio by region
child_ica %>% 
  ggplot(aes(x = factor(DID), y = factor(H002), fill = factor(H002))) +
  geom_bar(aes(x = factor(DID), y = factor(H002)), stat = "identity", position = "fill") +
  facet_grid(.~RID, scales = "free") +
  facet_wrap(.~RID, scales = "free") +
  theme(axis.text.x = element_text(angle = 90, size = 7, hjust = 1)) +
  coord_cartesian(ylim = c(0,1), clip = "off", expand = FALSE)

Gilgit-Baltistan. house type
child_ica %>% 
  filter(C003 != "NA", RID == 6) %>% 
  group_by(DID, H002) %>% 
  mutate(avg_enrollment = mean(C003 == 3)) %>% 
  ggplot(aes(factor(DID), avg_enrollment, group = H002, fill = factor(H002))) +
  geom_col(position = "dodge")

Hunza. house type
child_ica %>% left_join(house, by = "HHID") %>% 
  filter(DID.x == 266) %>% 
  ggplot(aes(H002.x)) +
  geom_histogram(aes(H002.x), stat = "count")
## Warning: Ignoring unknown parameters: binwidth, bins, pad

Average Number of Children in a household

child_ica %>% 
  group_by(HHID) %>% 
  mutate(n_children = length(unique(CID))) %>% 
  ungroup() %>% 
  group_by(DID) %>% 
  mutate(avg_n_children = mean(n_children)) %>% 
  ggplot(aes(DID, avg_n_children, color = factor(RID))) +
  geom_point()

child_ica %>% 
  ggplot(aes(C005)) +
  geom_histogram(stat = "count") +
  coord_flip()
## Warning: Ignoring unknown parameters: binwidth, bins, pad

Hunza, house type
child_ica %>% 
  filter(DID == 266) %>% 
  group_by(H002) %>% 
  summarize(n = n()) %>% 
  ggplot(aes(H002, n)) +
  geom_col() +
  geom_text(aes(label = n), vjust = 1.5, color = "white") +
  xlab("House Type") +
  ylab("Number of Households")
## `summarise()` ungrouping output (override with `.groups` argument)

Hunza, number of children in households
child_ica%>% 
  filter(DID == 266) %>% 
  group_by(HHID) %>%
  summarize(n = length(unique(CID))) %>% 
  ungroup() %>% 
  group_by(n) %>% 
  summarize(num = n()) %>% 
  ggplot(aes(factor(n), num)) +
  geom_col() +
  geom_text(aes(label = num), vjust = -0.5) +
  xlab("Number of Children in Each Household") +
  ylab("Number of Households") 
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

VID

Gilgit (DID == 260)
child_ica %>%
  filter(!is.na(C002), DID == 260) %>% 
  group_by(VID, C002) %>% 
  summarize(avg_enroll = mean(C003 == 3)) %>% 
  ggplot(aes(VID, avg_enroll, color = factor(C002))) +
  geom_point()
## `summarise()` regrouping output by 'VID' (override with `.groups` argument)

Preparation for Logistic Regression Analysis

Making Dataframe With Dummy Variables

child_ica_dummy <- child_ica %>% filter(!is.na(C002), !is.na(C003), !is.na(PR004), !is.na(PR009))
child_ica_dummy$C002_01 <- ifelse(child_ica_dummy$C002 == -1, 1, 0)
child_ica_dummy$C003_01 <- ifelse(child_ica_dummy$C003 == 3, 1, 0)
child_ica_dummy$PR004_01 <- ifelse(child_ica_dummy$PR004 == -1, 1, 0)
child_ica_dummy$PR009_01 <- ifelse(child_ica_dummy$PR009 == -1, 1, 0)
child_ica_dummy$PR004_PR009_01 <- ifelse(child_ica_dummy$PR004 == -1 | child_ica_dummy$PR009 == -1, 1, 0)
n_children_in_household <- child_ica_dummy %>% 
  group_by(HHID) %>% 
  summarize(n_children_in_household = length(unique(CID)))
## `summarise()` ungrouping output (override with `.groups` argument)
child_ica_dummy <- child_ica_dummy %>% left_join(n_children_in_household)
## Joining, by = "HHID"
child_ica_dummy$H002_1_01 <- ifelse(child_ica_dummy$H002 == 1, 1, 0)
child_ica_dummy$H002_2_01 <- ifelse(child_ica_dummy$H002 == 2, 1, 0)
child_ica_dummy$H002_3_01 <- ifelse(child_ica_dummy$H002 == 3, 1, 0)

Create dummy variable with regard to region

child_ica_dummy <- child_ica_dummy
child_ica_dummy$Panjab <- ifelse(child_ica_dummy$RID == 2, 1, 0)
child_ica_dummy$Sindh <- ifelse(child_ica_dummy$RID == 3, 1, 0)
child_ica_dummy$Balochistan <- ifelse(child_ica_dummy$RID == 4, 1, 0)
child_ica_dummy$Khyber_Pakhtunkhwa <- ifelse(child_ica_dummy$RID == 5, 1, 0)
child_ica_dummy$Gilgit_Baltistan <- ifelse(child_ica_dummy$RID == 6, 1, 0)
child_ica_dummy$Azad_Jammu_and_Kashmir <- ifelse(child_ica_dummy$RID == 7, 1, 0)
child_ica_dummy$Islamabad_ICT <- ifelse(child_ica_dummy$RID == 8, 1, 0)
child_ica_dummy$Federally_Administrated_Tribal_Areas <- ifelse(child_ica_dummy$RID == 9, 1, 0)

Eliminated NAs

child_ica %>% 
  summarize(C002_na = sum(is.na(C002)),
            C003_na = sum(is.na(C003)),
            PR004_na = sum(is.na(PR004)),
            PR009_na = sum(is.na(PR009)))

Eliminated Rows in Total

data.frame(original_rows = nrow(child_ica),
           eliminated_rows = nrow(child_ica) - nrow(child_ica_dummy),
           ratio = (nrow(child_ica)-nrow(child_ica_dummy))/nrow(child_ica))

Eliminated NAs Hunza

child_ica %>% 
  filter(DID == 266) %>% 
  summarize(C002_na = sum(is.na(C002)),
            C003_na = sum(is.na(C003)),
            PR004_na = sum(is.na(PR004)),
            PR009_na = sum(is.na(PR009)))

Eliminated Rows in Total Hunza

data.frame(original_rows = nrow(child_ica %>% filter(DID == 266)),
           eliminated_rows = nrow(child_ica %>% filter(DID == 266)) - nrow(child_ica_dummy %>% filter(DID == 266)),
           ratio = (nrow(child_ica %>% filter(DID == 266) %>% filter(DID == 266))-nrow(child_ica_dummy %>% filter(DID == 266)))/nrow(child_ica %>% filter(DID == 266)))

Generalized Linear Models

GLM C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01

glm_child <- glm(C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01, family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01, 
##     family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.8850  -1.2383   0.6563   0.8621   1.3071  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              0.934846   0.013780   67.84   <2e-16 ***
## n_children_in_household -0.055206   0.002935  -18.81   <2e-16 ***
## C002_01                 -0.572083   0.009041  -63.28   <2e-16 ***
## PR004_PR009_01           0.711645   0.009062   78.53   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 289072  on 245743  degrees of freedom
## AIC: 289080
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 490.13, df = 8, p-value < 2.2e-16

Hunza. C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01

glm_child_hunza <- glm(C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))

ggPredict(glm_child_hunza, se = TRUE, show.summary = TRUE, point = TRUE, colorAsFactor = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + C002_01 + PR004_PR009_01, 
##     family = "binomial", data = child_ica_dummy %>% filter(DID == 
##         266))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5255   0.3335   0.3611   0.4126   0.6056  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              2.74787    0.35406   7.761 8.43e-15 ***
## n_children_in_household -0.16354    0.07071  -2.313   0.0207 *  
## C002_01                  0.12227    0.19347   0.632   0.5274    
## PR004_PR009_01           0.44052    0.21451   2.054   0.0400 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 838.79  on 1609  degrees of freedom
## Residual deviance: 826.06  on 1606  degrees of freedom
## AIC: 834.06
## 
## Number of Fisher Scoring iterations: 5
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
Hosmer-Lemeshow
hoslem.test(x = glm_child_hunza$y, y = fitted(glm_child_hunza))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child_hunza$y, fitted(glm_child_hunza)
## X-squared = 9.9762, df = 8, p-value = 0.2667
exponential transformation
exp(glm_child_hunza$coefficients)
##             (Intercept) n_children_in_household                 C002_01 
##              15.6093953               0.8491308               1.1300567 
##          PR004_PR009_01 
##               1.5535098
confidence interval (intercept and coefficient)
confint(glm_child_hunza, level = 0.95)
## Waiting for profiling to be done...
##                              2.5 %      97.5 %
## (Intercept)              2.0690843  3.45821705
## n_children_in_household -0.3018419 -0.02433219
## C002_01                 -0.2575389  0.50262655
## PR004_PR009_01           0.0115443  0.85445608
exponential transformation of confidence interval (odds ratio of intercept and coefficient)
exp(confint(glm_child_hunza, level = 0.95))
## Waiting for profiling to be done...
##                             2.5 %     97.5 %
## (Intercept)             7.9175694 31.7602989
## n_children_in_household 0.7394550  0.9759614
## C002_01                 0.7729515  1.6530574
## PR004_PR009_01          1.0116112  2.3500958
AIC
extractAIC(glm_child_hunza)
## [1]   4.0000 834.0574
BIC
extractAIC(glm_child_hunza, k = log(nrow(glm_child_hunza$data)))
## [1]   4.0000 855.5934
effectiveness of explanatory variables
glm_child_hunza_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))
anova <- anova(glm_child_hunza_null, glm_child_hunza, test = "Chisq")
anova
variables selection
step(glm_child_hunza_null, direction = "both", 
     scope = (~ C001 + C002_01 + PR004_01 + PR009_01))
## Start:  AIC=840.79
## C003_01 ~ 1
## 
##            Df Deviance    AIC
## + C001      1   732.49 736.49
## + PR004_01  1   824.25 828.25
## + PR009_01  1   830.82 834.82
## <none>          838.79 840.79
## + C002_01   1   838.58 842.58
## 
## Step:  AIC=736.49
## C003_01 ~ C001
## 
##            Df Deviance    AIC
## + PR004_01  1   707.11 713.11
## + PR009_01  1   717.26 723.26
## <none>          732.49 736.49
## + C002_01   1   732.49 738.49
## - C001      1   838.79 840.79
## 
## Step:  AIC=713.11
## C003_01 ~ C001 + PR004_01
## 
##            Df Deviance    AIC
## <none>          707.11 713.11
## + PR009_01  1   705.74 713.74
## + C002_01   1   707.11 715.11
## - PR004_01  1   732.49 736.49
## - C001      1   824.25 828.25
## 
## Call:  glm(formula = C003_01 ~ C001 + PR004_01, family = "binomial", 
##     data = child_ica_dummy %>% filter(DID == 266))
## 
## Coefficients:
## (Intercept)         C001     PR004_01  
##     -0.4789       0.3039       1.0316  
## 
## Degrees of Freedom: 1609 Total (i.e. Null);  1607 Residual
## Null Deviance:       838.8 
## Residual Deviance: 707.1     AIC: 713.1
multicolinearity
vif(glm_child_hunza)
## n_children_in_household                 C002_01          PR004_PR009_01 
##                1.087196                1.006453                1.080729

GLM C003_01 ~ n_children_in_household + H002 + PR004_PR009_01

glm_child <- glm(C003_01 ~ n_children_in_household + H002 + PR004_PR009_01, family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + H002 + PR004_PR009_01, 
##     family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##    Min      1Q  Median      3Q     Max  
## -1.927  -1.297   0.723   0.869   1.186  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              0.094686   0.016614   5.699  1.2e-08 ***
## n_children_in_household -0.049487   0.002933 -16.871  < 2e-16 ***
## H002                     0.379659   0.006366  59.636  < 2e-16 ***
## PR004_PR009_01           0.502871   0.009477  53.061  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 289446  on 245743  degrees of freedom
## AIC: 289454
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 355.01, df = 8, p-value < 2.2e-16

Hunza. C003_01 ~ n_children_in_household + H002 + PR004_PR009_01

glm_child_hunza <- glm(C003_01 ~ n_children_in_household + H002 + PR004_PR009_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))

ggPredict(glm_child_hunza, se = TRUE, show.summary = TRUE, point = TRUE, colorAsFactor = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + H002 + PR004_PR009_01, 
##     family = "binomial", data = child_ica_dummy %>% filter(DID == 
##         266))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.6103   0.3376   0.3660   0.4038   0.5903  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              3.23266    0.49678   6.507 7.65e-11 ***
## n_children_in_household -0.16663    0.07002  -2.380   0.0173 *  
## H002                    -0.20362    0.16445  -1.238   0.2157    
## PR004_PR009_01           0.51067    0.22107   2.310   0.0209 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 838.79  on 1609  degrees of freedom
## Residual deviance: 824.91  on 1606  degrees of freedom
## AIC: 832.91
## 
## Number of Fisher Scoring iterations: 5
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
Hosmer-Lemeshow
hoslem.test(x = glm_child_hunza$y, y = fitted(glm_child_hunza))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child_hunza$y, fitted(glm_child_hunza)
## X-squared = 15.469, df = 8, p-value = 0.05064
exponential transformation
exp(glm_child_hunza$coefficients)
##             (Intercept) n_children_in_household                    H002 
##              25.3469307               0.8465115               0.8157746 
##          PR004_PR009_01 
##               1.6664083
confidence interval (intercept and coefficient)
confint(glm_child_hunza, level = 0.95)
## Waiting for profiling to be done...
##                               2.5 %      97.5 %
## (Intercept)              2.27327932  4.22193092
## n_children_in_household -0.30339827 -0.02854365
## H002                    -0.52934555  0.11610748
## PR004_PR009_01           0.06943756  0.93802587
exponential transformation of confidence interval (odds ratio of intercept and coefficient)
exp(confint(glm_child_hunza, level = 0.95))
## Waiting for profiling to be done...
##                             2.5 %     97.5 %
## (Intercept)             9.7111948 68.1649785
## n_children_in_household 0.7383050  0.9718599
## H002                    0.5889903  1.1231166
## PR004_PR009_01          1.0719051  2.5549327
AIC
extractAIC(glm_child_hunza)
## [1]   4.0000 832.9071
BIC
extractAIC(glm_child_hunza, k = log(nrow(glm_child_hunza$data)))
## [1]   4.0000 854.4431
effectiveness of explanatory variables
glm_child_hunza_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))
anova <- anova(glm_child_hunza_null, glm_child_hunza, test = "Chisq")
anova
variables selection
step(glm_child_hunza_null, direction = "both", 
     scope = (~ C001 + C002_01 + PR004_01 + PR009_01))
## Start:  AIC=840.79
## C003_01 ~ 1
## 
##            Df Deviance    AIC
## + C001      1   732.49 736.49
## + PR004_01  1   824.25 828.25
## + PR009_01  1   830.82 834.82
## <none>          838.79 840.79
## + C002_01   1   838.58 842.58
## 
## Step:  AIC=736.49
## C003_01 ~ C001
## 
##            Df Deviance    AIC
## + PR004_01  1   707.11 713.11
## + PR009_01  1   717.26 723.26
## <none>          732.49 736.49
## + C002_01   1   732.49 738.49
## - C001      1   838.79 840.79
## 
## Step:  AIC=713.11
## C003_01 ~ C001 + PR004_01
## 
##            Df Deviance    AIC
## <none>          707.11 713.11
## + PR009_01  1   705.74 713.74
## + C002_01   1   707.11 715.11
## - PR004_01  1   732.49 736.49
## - C001      1   824.25 828.25
## 
## Call:  glm(formula = C003_01 ~ C001 + PR004_01, family = "binomial", 
##     data = child_ica_dummy %>% filter(DID == 266))
## 
## Coefficients:
## (Intercept)         C001     PR004_01  
##     -0.4789       0.3039       1.0316  
## 
## Degrees of Freedom: 1609 Total (i.e. Null);  1607 Residual
## Null Deviance:       838.8 
## Residual Deviance: 707.1     AIC: 713.1
# step(glm_child_null, direction = "both", 
#      scope = (~ C001 + C002 + PR004_01 + PR009_01 + 
#                 Panjab + Sindh + Balochistan + Khyber_Pakhtunkhwa + Gilgit_Baltistan + 
#                 Azad_Jammu_and_Kashmir + Islamabad_ICT))
multicolinearity
vif(glm_child_hunza)
## n_children_in_household                    H002          PR004_PR009_01 
##                1.083416                1.088995                1.147110

GLM C003_01 ~ n_children_in_household + PR004_PR009_01 + H002

glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + H002, family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     H002, family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##    Min      1Q  Median      3Q     Max  
## -1.927  -1.297   0.723   0.869   1.186  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              0.094686   0.016614   5.699  1.2e-08 ***
## n_children_in_household -0.049487   0.002933 -16.871  < 2e-16 ***
## PR004_PR009_01           0.502871   0.009477  53.061  < 2e-16 ***
## H002                     0.379659   0.006366  59.636  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 289446  on 245743  degrees of freedom
## AIC: 289454
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 355.01, df = 8, p-value < 2.2e-16

Hunza. C003_01 ~ n_children_in_household + PR004_PR009_01 + H002

glm_child_hunza <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + H002, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))

ggPredict(glm_child_hunza, se = TRUE, show.summary = TRUE, point = TRUE, colorAsFactor = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     H002, family = "binomial", data = child_ica_dummy %>% filter(DID == 
##     266))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.6103   0.3376   0.3660   0.4038   0.5903  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              3.23266    0.49678   6.507 7.65e-11 ***
## n_children_in_household -0.16663    0.07002  -2.380   0.0173 *  
## PR004_PR009_01           0.51067    0.22107   2.310   0.0209 *  
## H002                    -0.20362    0.16445  -1.238   0.2157    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 838.79  on 1609  degrees of freedom
## Residual deviance: 824.91  on 1606  degrees of freedom
## AIC: 832.91
## 
## Number of Fisher Scoring iterations: 5
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 355.01, df = 8, p-value < 2.2e-16

Multicolinearity house type and parent edu?

child_ica_dummy %>% 
  group_by(H002) %>% 
  mutate(n_parents = length(PR004_PR009_01))%>% 
  ggplot(aes(H002, n_parents, fill = factor(PR004_PR009_01))) +
  geom_col(position = "stack") +
  xlab("Type of Houses") +
  ylab("Number of Parents") +
  labs(fill = "Parents Education\n(at least one of them, primary or more)")
The more educated the parents are, the better house they are living in.

GLM C003_01 ~ n_children_in_household + H002 + C002_01

glm_child <- glm(C003_01 ~ n_children_in_household + H002 + C002_01, family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + H002 + C002_01, 
##     family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.0122  -1.2391   0.7184   0.8909   1.2724  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              0.424111   0.016894   25.10   <2e-16 ***
## n_children_in_household -0.052085   0.002941  -17.71   <2e-16 ***
## H002                     0.503563   0.006129   82.16   <2e-16 ***
## C002_01                 -0.575354   0.009055  -63.54   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 288199  on 245743  degrees of freedom
## AIC: 288207
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 818.15, df = 8, p-value < 2.2e-16

Hunza. C003_01 ~ n_children_in_household + H002 + C002_01

glm_child_hunza <- glm(C003_01 ~ n_children_in_household + H002 + C002_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))

ggPredict(glm_child_hunza, se = TRUE, show.summary = TRUE, point = TRUE, colorAsFactor = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + H002 + C002_01, 
##     family = "binomial", data = child_ica_dummy %>% filter(DID == 
##         266))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5804   0.3373   0.3731   0.4124   0.5843  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              3.48657    0.48393   7.205 5.82e-13 ***
## n_children_in_household -0.20822    0.06887  -3.023   0.0025 ** 
## H002                    -0.11589    0.15847  -0.731   0.4646    
## C002_01                  0.13022    0.19320   0.674   0.5003    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 838.79  on 1609  degrees of freedom
## Residual deviance: 829.57  on 1606  degrees of freedom
## AIC: 837.57
## 
## Number of Fisher Scoring iterations: 5
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child_hunza$y, y = fitted(glm_child_hunza))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child_hunza$y, fitted(glm_child_hunza)
## X-squared = 12.015, df = 8, p-value = 0.1505
exponential transformation
exp(glm_child_hunza$coefficients)
##             (Intercept) n_children_in_household                    H002 
##              32.6735833               0.8120302               0.8905720 
##                 C002_01 
##               1.1390753
confidence interval (intercept and coefficient)
confint(glm_child_hunza, level = 0.95)
## Waiting for profiling to be done...
##                              2.5 %      97.5 %
## (Intercept)              2.5525708  4.45124970
## n_children_in_household -0.3424699 -0.07212852
## H002                    -0.4291928  0.19275626
## C002_01                 -0.2490334  0.51006070
exponential transformation of confidence interval (odds ratio of intercept and coefficient)
exp(confint(glm_child_hunza, level = 0.95))
## Waiting for profiling to be done...
##                              2.5 %     97.5 %
## (Intercept)             12.8400705 85.7340185
## n_children_in_household  0.7100145  0.9304113
## H002                     0.6510344  1.2125872
## C002_01                  0.7795540  1.6653923
AIC
extractAIC(glm_child_hunza)
## [1]   4.0000 837.5663
BIC
extractAIC(glm_child_hunza, k = log(nrow(glm_child_hunza$data)))
## [1]   4.0000 859.1022
effectiveness of explanatory variables
glm_child_hunza_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))
anova <- anova(glm_child_hunza_null, glm_child_hunza, test = "Chisq")
anova
variables selection
step(glm_child_hunza_null, direction = "both", 
     scope = (~ C001 + C002_01 + PR004_01 + PR009_01))
## Start:  AIC=840.79
## C003_01 ~ 1
## 
##            Df Deviance    AIC
## + C001      1   732.49 736.49
## + PR004_01  1   824.25 828.25
## + PR009_01  1   830.82 834.82
## <none>          838.79 840.79
## + C002_01   1   838.58 842.58
## 
## Step:  AIC=736.49
## C003_01 ~ C001
## 
##            Df Deviance    AIC
## + PR004_01  1   707.11 713.11
## + PR009_01  1   717.26 723.26
## <none>          732.49 736.49
## + C002_01   1   732.49 738.49
## - C001      1   838.79 840.79
## 
## Step:  AIC=713.11
## C003_01 ~ C001 + PR004_01
## 
##            Df Deviance    AIC
## <none>          707.11 713.11
## + PR009_01  1   705.74 713.74
## + C002_01   1   707.11 715.11
## - PR004_01  1   732.49 736.49
## - C001      1   824.25 828.25
## 
## Call:  glm(formula = C003_01 ~ C001 + PR004_01, family = "binomial", 
##     data = child_ica_dummy %>% filter(DID == 266))
## 
## Coefficients:
## (Intercept)         C001     PR004_01  
##     -0.4789       0.3039       1.0316  
## 
## Degrees of Freedom: 1609 Total (i.e. Null);  1607 Residual
## Null Deviance:       838.8 
## Residual Deviance: 707.1     AIC: 713.1
# step(glm_child_null, direction = "both", 
#      scope = (~ C001 + C002 + PR004_01 + PR009_01 + 
#                 Panjab + Sindh + Balochistan + Khyber_Pakhtunkhwa + Gilgit_Baltistan + 
#                 Azad_Jammu_and_Kashmir + Islamabad_ICT))
multicolinearity
vif(glm_child_hunza)
## n_children_in_household                    H002                 C002_01 
##                1.021627                1.015482                1.006245

GLM C003_01 ~ n_children_in_household + C002_01 + H002

glm_child <- glm(C003_01 ~ n_children_in_household + C002_01 + H002, family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + C002_01 + H002, 
##     family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.0122  -1.2391   0.7184   0.8909   1.2724  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              0.424111   0.016894   25.10   <2e-16 ***
## n_children_in_household -0.052085   0.002941  -17.71   <2e-16 ***
## C002_01                 -0.575354   0.009055  -63.54   <2e-16 ***
## H002                     0.503563   0.006129   82.16   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 288199  on 245743  degrees of freedom
## AIC: 288207
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 818.15, df = 8, p-value < 2.2e-16

Hunza. C003_01 ~n_children_in_household + C002_01 + H002

glm_child_hunza <- glm(C003_01 ~n_children_in_household + C002_01 + H002, family = "binomial", data = child_ica_dummy %>% filter(DID == 266))

ggPredict(glm_child_hunza, se = TRUE, show.summary = TRUE, point = TRUE, colorAsFactor = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + C002_01 + H002, 
##     family = "binomial", data = child_ica_dummy %>% filter(DID == 
##         266))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5804   0.3373   0.3731   0.4124   0.5843  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              3.48657    0.48393   7.205 5.82e-13 ***
## n_children_in_household -0.20822    0.06887  -3.023   0.0025 ** 
## C002_01                  0.13022    0.19320   0.674   0.5003    
## H002                    -0.11589    0.15847  -0.731   0.4646    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 838.79  on 1609  degrees of freedom
## Residual deviance: 829.57  on 1606  degrees of freedom
## AIC: 837.57
## 
## Number of Fisher Scoring iterations: 5
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child_hunza$y, y = fitted(glm_child_hunza))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child_hunza$y, fitted(glm_child_hunza)
## X-squared = 12.015, df = 8, p-value = 0.1505

C002_01 chi Districts in Gilgit-Baltistan

child_ica_dummy %>% filter(RID == 6) %>% summarize(unique(DID))
sapply(260:266, function(dist){
  glm_child_dist <- glm(C003_01 ~ C002_01, family = "binomial", data = child_ica_dummy %>% filter(RID == 6, DID == dist))
  glm_child_dist_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy %>% filter(RID == 6, DID == dist))
anova <- anova(glm_child_dist_null, glm_child_dist, test = "Chisq")
data.frame(result = anova$`Pr(>Chi)`[2], DID = dist, chi_0.05 = anova$`Pr(>Chi)`[2] >= 0.05)
})
##          [,1]       [,2]          [,3]       [,4]       [,5]      [,6]     
## result   0.02323511 1.711351e-125 0.02518717 0.04226537 0.1507538 0.7348869
## DID      260        261           262        263        264       265      
## chi_0.05 FALSE      FALSE         FALSE      FALSE      TRUE      TRUE     
##          [,7]     
## result   0.6463022
## DID      266      
## chi_0.05 TRUE
Within Gilgit-Baltistan, 3 of 7 districts including Hunza has Chi > 0.05

GLM C003_01 ~ C002_01 and Anova for All Districts

DID_unique <- child_ica_dummy %>% summarize(unique(DID))
DID_unique <- as.matrix(DID_unique)
DID_unique <- as.numeric(DID_unique[,1])

dist_logist <- sapply(DID_unique, function(dist){
  
  glm_child_dist <- glm(C003_01 ~ C002_01, 
                        family = "binomial", 
                        data = child_ica_dummy %>% 
                          filter(DID == dist))
  
  glm_child_dist_null <- glm(C003_01~1, 
                             family = "binomial", 
                             data = child_ica_dummy %>% 
                               filter(DID == dist))

  anova <- anova(glm_child_dist_null, glm_child_dist, test = "Chisq")

data.frame(DID = dist, 
           chi = anova$`Pr(>Chi)`[2], 
           chi_largerThan_0.05 = anova$`Pr(>Chi)`[2] >= 0.05, 
           chi_largerThan_0.1 = anova$`Pr(>Chi)`[2] >= 0.1)
})

anova_dist <- as.data.frame(as.tibble(t(dist_logist)))
## Warning: `as.tibble()` is deprecated as of tibble 2.0.0.
## Please use `as_tibble()` instead.
## The signature and semantics have changed, see `?as_tibble`.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.
anova_dist$DID <- as.numeric(anova_dist$DID)
anova_dist$chi <- as.numeric(anova_dist$chi)
anova_dist$chi_largerThan_0.1 <- as.logical(anova_dist$chi_largerThan_0.1)
anova_dist$chi_largerThan_0.05 <- as.logical(anova_dist$chi_largerThan_0.05)
anova_dist %>% 
  summarize(total_chi_largerThan_0.05 = sum(as.numeric(chi_largerThan_0.05)),
            ratio_0.05 = sum(as.numeric(chi_largerThan_0.05))/length(chi_largerThan_0.05), 
            total_chi_largerThan_0.1 = sum(as.numeric(chi_largerThan_0.1)),
            ratio_0.1 = sum(as.numeric(chi_largerThan_0.1))/length(chi_largerThan_0.1))
DIDs with chi >= 0.05
anova_dist %>% 
  filter(chi_largerThan_0.05 == TRUE) %>% 
  select(-chi_largerThan_0.1)
DID chi > 0.05
DID_chi_0.05 <- anova_dist %>% 
  filter(chi_largerThan_0.05 == TRUE) %>% 
  select(DID) %>% 
  pull()
DID_chi_0.05
##  [1] 146 151 156 158 159 162 163 164 167 171 173 176 178 179 189 196 245 257 264
## [20] 265 266 267 268 269 270 271 272 273 274 276
DIDs with chi >= 0.1
anova_dist %>% 
  filter(chi_largerThan_0.1 == TRUE) %>% 
  select(-chi_largerThan_0.05)
DID chi > 0.1
DID_chi_0.1 <- anova_dist %>% 
  filter(chi_largerThan_0.1 == TRUE) %>% 
  select(DID) %>% 
  pull()
DID with C002 == -1, chi > 0.05
child %>% 
  filter(C002 == -1) %>% 
  group_by(DID) %>% 
  mutate(avg_learning = mean(C010, na.rm = TRUE)) %>% 
  ggplot(aes(DID, avg_learning)) +
  geom_point(data = child %>% 
               filter(C002 == -1 & DID == DID_chi_0.05) %>% 
               group_by(DID) %>%
               mutate(avg_learning = mean(C010, na.rm = TRUE)), size = 4, color = "red") +
  geom_point(aes(color = RID)) +
  geom_text(aes(label = DID), nudge_x = 5, check_overlap = TRUE)
## Warning in DID == DID_chi_0.05: 長いオブジェクトの長さが短いオブジェクトの長さの
## 倍数になっていません

C001 & enroll ; DID chi > 0.1
child %>% 
  filter(DID == DID_chi_0.1) %>%
  group_by(DID, C001) %>% 
  mutate(enroll = mean(C003 == 3)) %>% 
  ggplot(aes(C001, enroll)) +
  geom_point() +
  geom_smooth() +
  facet_wrap(DID~.)
## Warning in DID == DID_chi_0.1: 長いオブジェクトの長さが短いオブジェクトの長さの
## 倍数になっていません
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

GLM C003_01 ~ n_children_in_household + H002 + C002_01 + as.factor(DID)

glm_child <- glm(C003_01 ~ n_children_in_household + H002 + C002_01 + as.factor(DID), family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + H002 + C002_01 + 
##     as.factor(DID), family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5433  -1.0857   0.6008   0.8455   1.6963  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              1.856387   0.090069  20.611  < 2e-16 ***
## n_children_in_household -0.030048   0.003248  -9.252  < 2e-16 ***
## H002                     0.238646   0.007614  31.342  < 2e-16 ***
## C002_01                 -0.661321   0.009514 -69.509  < 2e-16 ***
## as.factor(DID)147       -0.531421   0.106816  -4.975 6.52e-07 ***
## as.factor(DID)148       -0.080366   0.111073  -0.724 0.469342    
## as.factor(DID)149       -0.729406   0.102725  -7.101 1.24e-12 ***
## as.factor(DID)150       -0.755588   0.101817  -7.421 1.16e-13 ***
## as.factor(DID)151       -0.188233   0.111184  -1.693 0.090458 .  
## as.factor(DID)152       -0.294774   0.115283  -2.557 0.010559 *  
## as.factor(DID)153       -0.876057   0.103921  -8.430  < 2e-16 ***
## as.factor(DID)154       -0.635280   0.105517  -6.021 1.74e-09 ***
## as.factor(DID)155       -0.824648   0.101025  -8.163 3.27e-16 ***
## as.factor(DID)156       -1.295679   0.115699 -11.199  < 2e-16 ***
## as.factor(DID)157       -0.457369   0.110809  -4.128 3.67e-05 ***
## as.factor(DID)158       -0.353198   0.107855  -3.275 0.001058 ** 
## as.factor(DID)159       -1.312255   0.106407 -12.332  < 2e-16 ***
## as.factor(DID)160       -0.699375   0.104766  -6.676 2.46e-11 ***
## as.factor(DID)161       -1.308519   0.097458 -13.426  < 2e-16 ***
## as.factor(DID)162        0.210991   0.125085   1.687 0.091645 .  
## as.factor(DID)163        0.071172   0.127802   0.557 0.577599    
## as.factor(DID)164       -0.508366   0.113661  -4.473 7.73e-06 ***
## as.factor(DID)165       -0.991673   0.104159  -9.521  < 2e-16 ***
## as.factor(DID)166       -0.869402   0.104533  -8.317  < 2e-16 ***
## as.factor(DID)167        0.526340   0.130238   4.041 5.31e-05 ***
## as.factor(DID)169       -1.655333   0.099225 -16.683  < 2e-16 ***
## as.factor(DID)170       -0.488745   0.108583  -4.501 6.76e-06 ***
## as.factor(DID)171       -0.436012   0.110567  -3.943 8.03e-05 ***
## as.factor(DID)172       -0.624901   0.106396  -5.873 4.27e-09 ***
## as.factor(DID)173       -0.306239   0.109770  -2.790 0.005274 ** 
## as.factor(DID)174       -0.895481   0.105491  -8.489  < 2e-16 ***
## as.factor(DID)175       -1.320879   0.101365 -13.031  < 2e-16 ***
## as.factor(DID)176       -0.464176   0.114670  -4.048 5.17e-05 ***
## as.factor(DID)177       -0.853066   0.102982  -8.284  < 2e-16 ***
## as.factor(DID)178        0.182187   0.117181   1.555 0.120006    
## as.factor(DID)179       -0.718075   0.105393  -6.813 9.54e-12 ***
## as.factor(DID)180       -0.673160   0.105975  -6.352 2.12e-10 ***
## as.factor(DID)181       -0.158811   0.115371  -1.377 0.168658    
## as.factor(DID)182       -1.252132   0.100397 -12.472  < 2e-16 ***
## as.factor(DID)183       -0.624218   0.107630  -5.800 6.64e-09 ***
## as.factor(DID)184       -1.279621   0.099900 -12.809  < 2e-16 ***
## as.factor(DID)185        0.315727   0.116482   2.711 0.006718 ** 
## as.factor(DID)186       -1.527947   0.099616 -15.338  < 2e-16 ***
## as.factor(DID)187       -1.365833   0.101106 -13.509  < 2e-16 ***
## as.factor(DID)188       -1.201238   0.100289 -11.978  < 2e-16 ***
## as.factor(DID)189       -0.279213   0.108191  -2.581 0.009858 ** 
## as.factor(DID)190       -0.984810   0.108473  -9.079  < 2e-16 ***
## as.factor(DID)191       -0.856359   0.108294  -7.908 2.62e-15 ***
## as.factor(DID)192       -0.907494   0.102627  -8.843  < 2e-16 ***
## as.factor(DID)193       -0.062637   0.110415  -0.567 0.570519    
## as.factor(DID)194       -1.641845   0.099810 -16.450  < 2e-16 ***
## as.factor(DID)195       -1.800501   0.100992 -17.828  < 2e-16 ***
## as.factor(DID)196       -0.643909   0.109777  -5.866 4.47e-09 ***
## as.factor(DID)197       -1.333062   0.100752 -13.231  < 2e-16 ***
## as.factor(DID)198       -1.750612   0.099615 -17.574  < 2e-16 ***
## as.factor(DID)199        0.411202   0.116763   3.522 0.000429 ***
## as.factor(DID)200       -1.504418   0.102715 -14.647  < 2e-16 ***
## as.factor(DID)202       -1.430222   0.100851 -14.182  < 2e-16 ***
## as.factor(DID)203       -1.244656   0.099211 -12.546  < 2e-16 ***
## as.factor(DID)204       -0.864444   0.103348  -8.364  < 2e-16 ***
## as.factor(DID)205       -1.475565   0.098277 -15.014  < 2e-16 ***
## as.factor(DID)206       -1.502885   0.098520 -15.255  < 2e-16 ***
## as.factor(DID)207       -1.931182   0.100015 -19.309  < 2e-16 ***
## as.factor(DID)208       -1.433055   0.099404 -14.417  < 2e-16 ***
## as.factor(DID)209       -0.835136   0.108254  -7.715 1.21e-14 ***
## as.factor(DID)210       -1.656130   0.097887 -16.919  < 2e-16 ***
## as.factor(DID)211       -0.940851   0.107889  -8.721  < 2e-16 ***
## as.factor(DID)212       -1.184811   0.102333 -11.578  < 2e-16 ***
## as.factor(DID)213       -1.492092   0.098622 -15.129  < 2e-16 ***
## as.factor(DID)214       -1.980250   0.099212 -19.960  < 2e-16 ***
## as.factor(DID)215       -1.277002   0.096980 -13.168  < 2e-16 ***
## as.factor(DID)216       -1.475047   0.101753 -14.496  < 2e-16 ***
## as.factor(DID)217       -1.645835   0.099970 -16.463  < 2e-16 ***
## as.factor(DID)218       -1.391871   0.097592 -14.262  < 2e-16 ***
## as.factor(DID)219       -1.414559   0.100711 -14.046  < 2e-16 ***
## as.factor(DID)220       -2.067994   0.096336 -21.466  < 2e-16 ***
## as.factor(DID)221       -2.314209   0.101510 -22.798  < 2e-16 ***
## as.factor(DID)222       -1.569588   0.097845 -16.042  < 2e-16 ***
## as.factor(DID)223       -1.497175   0.105467 -14.196  < 2e-16 ***
## as.factor(DID)224       -2.024431   0.099120 -20.424  < 2e-16 ***
## as.factor(DID)225       -1.626735   0.099947 -16.276  < 2e-16 ***
## as.factor(DID)226       -1.380602   0.102822 -13.427  < 2e-16 ***
## as.factor(DID)227       -1.587866   0.098118 -16.183  < 2e-16 ***
## as.factor(DID)228       -2.149597   0.103253 -20.819  < 2e-16 ***
## as.factor(DID)229       -2.016092   0.101742 -19.816  < 2e-16 ***
## as.factor(DID)230       -1.952694   0.099633 -19.599  < 2e-16 ***
## as.factor(DID)231       -1.323551   0.098863 -13.388  < 2e-16 ***
## as.factor(DID)232       -0.149362   0.107653  -1.387 0.165307    
## as.factor(DID)233       -2.172729   0.104344 -20.823  < 2e-16 ***
## as.factor(DID)234       -1.481327   0.105267 -14.072  < 2e-16 ***
## as.factor(DID)235       -0.728578   0.102880  -7.082 1.42e-12 ***
## as.factor(DID)236       -1.011111   0.103059  -9.811  < 2e-16 ***
## as.factor(DID)237       -0.768246   0.105928  -7.253 4.09e-13 ***
## as.factor(DID)238        0.092540   0.116242   0.796 0.425973    
## as.factor(DID)239       -0.895697   0.102414  -8.746  < 2e-16 ***
## as.factor(DID)240       -1.139789   0.102174 -11.155  < 2e-16 ***
## as.factor(DID)241       -2.331117   0.098219 -23.734  < 2e-16 ***
## as.factor(DID)242       -1.090313   0.103367 -10.548  < 2e-16 ***
## as.factor(DID)243       -1.092509   0.108519 -10.067  < 2e-16 ***
## as.factor(DID)244       -0.389501   0.110870  -3.513 0.000443 ***
## as.factor(DID)245        0.307911   0.125277   2.458 0.013978 *  
## as.factor(DID)246       -1.693232   0.100045 -16.925  < 2e-16 ***
## as.factor(DID)247       -1.233913   0.107645 -11.463  < 2e-16 ***
## as.factor(DID)248       -0.789366   0.100028  -7.891 2.99e-15 ***
## as.factor(DID)249       -0.023336   0.124757  -0.187 0.851620    
## as.factor(DID)250       -0.687107   0.105155  -6.534 6.39e-11 ***
## as.factor(DID)251       -0.540857   0.106874  -5.061 4.18e-07 ***
## as.factor(DID)252       -0.860171   0.103757  -8.290  < 2e-16 ***
## as.factor(DID)253       -0.422480   0.111991  -3.772 0.000162 ***
## as.factor(DID)254       -1.716154   0.115165 -14.902  < 2e-16 ***
## as.factor(DID)255       -0.718627   0.108954  -6.596 4.23e-11 ***
## as.factor(DID)256       -0.192197   0.106990  -1.796 0.072430 .  
## as.factor(DID)257        0.525470   0.140459   3.741 0.000183 ***
## as.factor(DID)258       -1.393632   0.099084 -14.065  < 2e-16 ***
## as.factor(DID)259       -0.803014   0.107179  -7.492 6.77e-14 ***
## as.factor(DID)260       -0.030875   0.108454  -0.285 0.775887    
## as.factor(DID)261       -2.203777   0.098347 -22.408  < 2e-16 ***
## as.factor(DID)262       -0.562816   0.103736  -5.425 5.78e-08 ***
## as.factor(DID)263       -0.168283   0.107682  -1.563 0.118105    
## as.factor(DID)264       -0.378301   0.105640  -3.581 0.000342 ***
## as.factor(DID)265        0.045178   0.114401   0.395 0.692909    
## as.factor(DID)266        0.651750   0.130083   5.010 5.44e-07 ***
## as.factor(DID)267       -0.138255   0.106774  -1.295 0.195374    
## as.factor(DID)268       -0.998930   0.104179  -9.589  < 2e-16 ***
## as.factor(DID)269        0.089197   0.112064   0.796 0.426058    
## as.factor(DID)270        0.047179   0.111679   0.422 0.672695    
## as.factor(DID)271        0.316533   0.115617   2.738 0.006186 ** 
## as.factor(DID)272        0.145384   0.112480   1.293 0.196175    
## as.factor(DID)273       -0.121213   0.109283  -1.109 0.267356    
## as.factor(DID)274        0.006667   0.120782   0.055 0.955983    
## as.factor(DID)275       -1.031188   0.101880 -10.122  < 2e-16 ***
## as.factor(DID)276       -0.061414   0.106435  -0.577 0.563931    
## as.factor(DID)277       -0.253596   0.146904  -1.726 0.084300 .  
## as.factor(DID)278       -0.665286   0.098548  -6.751 1.47e-11 ***
## as.factor(DID)279       -1.754758   0.100067 -17.536  < 2e-16 ***
## as.factor(DID)280       -0.906296   0.102311  -8.858  < 2e-16 ***
## as.factor(DID)281       -0.667260   0.101980  -6.543 6.03e-11 ***
## as.factor(DID)282       -0.588057   0.109290  -5.381 7.42e-08 ***
## as.factor(DID)284       -0.676702   0.104745  -6.460 1.04e-10 ***
## as.factor(DID)287       -0.907293   0.104730  -8.663  < 2e-16 ***
## as.factor(DID)289       -0.917523   0.101097  -9.076  < 2e-16 ***
## as.factor(DID)290       -0.888790   0.104415  -8.512  < 2e-16 ***
## as.factor(DID)315       -0.736484   0.115424  -6.381 1.76e-10 ***
## as.factor(DID)316       -1.195621   0.103364 -11.567  < 2e-16 ***
## as.factor(DID)318       -1.593066   0.105466 -15.105  < 2e-16 ***
## as.factor(DID)319       -2.132163   0.099025 -21.532  < 2e-16 ***
## as.factor(DID)320       -1.398009   0.104773 -13.343  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 269755  on 245600  degrees of freedom
## AIC: 270049
## 
## Number of Fisher Scoring iterations: 5
## Warning in ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary =
## TRUE, : maximum three independent variables are allowed
## NULL
Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 330.28, df = 8, p-value < 2.2e-16

GLM C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(DID)

glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(DID), family = "binomial", data = child_ica_dummy)

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01 + as.factor(DID), family = "binomial", data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5059  -1.0668   0.6014   0.8352   1.7174  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              2.161512   0.088576  24.403  < 2e-16 ***
## n_children_in_household -0.028508   0.003251  -8.768  < 2e-16 ***
## PR004_PR009_01           0.347071   0.010500  33.054  < 2e-16 ***
## C002_01                 -0.664193   0.009519 -69.776  < 2e-16 ***
## as.factor(DID)147       -0.607598   0.106738  -5.692 1.25e-08 ***
## as.factor(DID)148       -0.093994   0.111048  -0.846 0.397314    
## as.factor(DID)149       -0.804572   0.102678  -7.836 4.65e-15 ***
## as.factor(DID)150       -0.784029   0.101828  -7.700 1.37e-14 ***
## as.factor(DID)151       -0.190780   0.111171  -1.716 0.086146 .  
## as.factor(DID)152       -0.334471   0.115245  -2.902 0.003705 ** 
## as.factor(DID)153       -0.900790   0.103870  -8.672  < 2e-16 ***
## as.factor(DID)154       -0.678256   0.105494  -6.429 1.28e-10 ***
## as.factor(DID)155       -0.842393   0.101011  -8.340  < 2e-16 ***
## as.factor(DID)156       -1.313467   0.115592 -11.363  < 2e-16 ***
## as.factor(DID)157       -0.452119   0.110799  -4.081 4.49e-05 ***
## as.factor(DID)158       -0.343909   0.107871  -3.188 0.001432 ** 
## as.factor(DID)159       -1.309531   0.106342 -12.314  < 2e-16 ***
## as.factor(DID)160       -0.738367   0.104702  -7.052 1.76e-12 ***
## as.factor(DID)161       -1.277624   0.097450 -13.111  < 2e-16 ***
## as.factor(DID)162        0.140673   0.125046   1.125 0.260600    
## as.factor(DID)163        0.076054   0.127784   0.595 0.551728    
## as.factor(DID)164       -0.566415   0.113648  -4.984 6.23e-07 ***
## as.factor(DID)165       -0.978953   0.104104  -9.404  < 2e-16 ***
## as.factor(DID)166       -0.830884   0.104540  -7.948 1.90e-15 ***
## as.factor(DID)167        0.509703   0.130242   3.914 9.10e-05 ***
## as.factor(DID)169       -1.723715   0.099078 -17.398  < 2e-16 ***
## as.factor(DID)170       -0.472852   0.108538  -4.357 1.32e-05 ***
## as.factor(DID)171       -0.437550   0.110541  -3.958 7.55e-05 ***
## as.factor(DID)172       -0.638525   0.106358  -6.004 1.93e-09 ***
## as.factor(DID)173       -0.267159   0.109787  -2.433 0.014957 *  
## as.factor(DID)174       -0.956716   0.105408  -9.076  < 2e-16 ***
## as.factor(DID)175       -1.331763   0.101311 -13.145  < 2e-16 ***
## as.factor(DID)176       -0.472402   0.114611  -4.122 3.76e-05 ***
## as.factor(DID)177       -0.852129   0.102954  -8.277  < 2e-16 ***
## as.factor(DID)178        0.202014   0.117141   1.725 0.084612 .  
## as.factor(DID)179       -0.718205   0.105390  -6.815 9.44e-12 ***
## as.factor(DID)180       -0.678592   0.105962  -6.404 1.51e-10 ***
## as.factor(DID)181       -0.117571   0.115344  -1.019 0.308054    
## as.factor(DID)182       -1.388829   0.100242 -13.855  < 2e-16 ***
## as.factor(DID)183       -0.659636   0.107584  -6.131 8.71e-10 ***
## as.factor(DID)184       -1.307346   0.099855 -13.092  < 2e-16 ***
## as.factor(DID)185        0.289557   0.116541   2.485 0.012970 *  
## as.factor(DID)186       -1.693746   0.099333 -17.051  < 2e-16 ***
## as.factor(DID)187       -1.457246   0.100967 -14.433  < 2e-16 ***
## as.factor(DID)188       -1.301699   0.100182 -12.993  < 2e-16 ***
## as.factor(DID)189       -0.357912   0.108160  -3.309 0.000936 ***
## as.factor(DID)190       -0.994167   0.108433  -9.168  < 2e-16 ***
## as.factor(DID)191       -0.905377   0.108243  -8.364  < 2e-16 ***
## as.factor(DID)192       -1.059675   0.102504 -10.338  < 2e-16 ***
## as.factor(DID)193       -0.060710   0.110471  -0.550 0.582622    
## as.factor(DID)194       -1.683985   0.099732 -16.885  < 2e-16 ***
## as.factor(DID)195       -1.882723   0.100866 -18.666  < 2e-16 ***
## as.factor(DID)196       -0.698907   0.109777  -6.367 1.93e-10 ***
## as.factor(DID)197       -1.429311   0.100649 -14.201  < 2e-16 ***
## as.factor(DID)198       -1.785280   0.099553 -17.933  < 2e-16 ***
## as.factor(DID)199        0.351822   0.116662   3.016 0.002564 ** 
## as.factor(DID)200       -1.457248   0.102727 -14.186  < 2e-16 ***
## as.factor(DID)202       -1.476437   0.100767 -14.652  < 2e-16 ***
## as.factor(DID)203       -1.442589   0.099023 -14.568  < 2e-16 ***
## as.factor(DID)204       -0.894943   0.103239  -8.669  < 2e-16 ***
## as.factor(DID)205       -1.587176   0.098221 -16.159  < 2e-16 ***
## as.factor(DID)206       -1.652157   0.098176 -16.828  < 2e-16 ***
## as.factor(DID)207       -2.029097   0.099744 -20.343  < 2e-16 ***
## as.factor(DID)208       -1.513908   0.099299 -15.246  < 2e-16 ***
## as.factor(DID)209       -1.027678   0.107955  -9.519  < 2e-16 ***
## as.factor(DID)210       -1.723325   0.097692 -17.640  < 2e-16 ***
## as.factor(DID)211       -1.029961   0.107737  -9.560  < 2e-16 ***
## as.factor(DID)212       -1.344713   0.102072 -13.174  < 2e-16 ***
## as.factor(DID)213       -1.582224   0.098485 -16.066  < 2e-16 ***
## as.factor(DID)214       -2.095527   0.099005 -21.166  < 2e-16 ***
## as.factor(DID)215       -1.359969   0.096766 -14.054  < 2e-16 ***
## as.factor(DID)216       -1.675878   0.101392 -16.529  < 2e-16 ***
## as.factor(DID)217       -1.709895   0.099764 -17.139  < 2e-16 ***
## as.factor(DID)218       -1.634273   0.097160 -16.820  < 2e-16 ***
## as.factor(DID)219       -1.542831   0.100494 -15.352  < 2e-16 ***
## as.factor(DID)220       -2.103022   0.096242 -21.851  < 2e-16 ***
## as.factor(DID)221       -2.455716   0.101216 -24.262  < 2e-16 ***
## as.factor(DID)222       -1.692441   0.097663 -17.329  < 2e-16 ***
## as.factor(DID)223       -1.604057   0.105314 -15.231  < 2e-16 ***
## as.factor(DID)224       -2.148013   0.098930 -21.712  < 2e-16 ***
## as.factor(DID)225       -1.852791   0.099548 -18.612  < 2e-16 ***
## as.factor(DID)226       -1.388128   0.102732 -13.512  < 2e-16 ***
## as.factor(DID)227       -1.677606   0.097855 -17.144  < 2e-16 ***
## as.factor(DID)228       -2.257037   0.102988 -21.916  < 2e-16 ***
## as.factor(DID)229       -1.971458   0.101725 -19.380  < 2e-16 ***
## as.factor(DID)230       -2.137518   0.099303 -21.525  < 2e-16 ***
## as.factor(DID)231       -1.403955   0.098699 -14.225  < 2e-16 ***
## as.factor(DID)232       -0.293014   0.107407  -2.728 0.006370 ** 
## as.factor(DID)233       -2.270591   0.104146 -21.802  < 2e-16 ***
## as.factor(DID)234       -1.583978   0.105021 -15.082  < 2e-16 ***
## as.factor(DID)235       -0.832266   0.102738  -8.101 5.46e-16 ***
## as.factor(DID)236       -1.105162   0.102937 -10.736  < 2e-16 ***
## as.factor(DID)237       -0.764733   0.105896  -7.222 5.14e-13 ***
## as.factor(DID)238        0.248673   0.116251   2.139 0.032427 *  
## as.factor(DID)239       -0.895892   0.102392  -8.750  < 2e-16 ***
## as.factor(DID)240       -1.242176   0.102128 -12.163  < 2e-16 ***
## as.factor(DID)241       -2.369040   0.098180 -24.129  < 2e-16 ***
## as.factor(DID)242       -1.101073   0.103317 -10.657  < 2e-16 ***
## as.factor(DID)243       -1.057710   0.108492  -9.749  < 2e-16 ***
## as.factor(DID)244       -0.560656   0.110773  -5.061 4.16e-07 ***
## as.factor(DID)245        0.261619   0.125199   2.090 0.036652 *  
## as.factor(DID)246       -1.706840   0.100013 -17.066  < 2e-16 ***
## as.factor(DID)247       -1.234545   0.107606 -11.473  < 2e-16 ***
## as.factor(DID)248       -0.945059   0.099728  -9.476  < 2e-16 ***
## as.factor(DID)249        0.026476   0.124731   0.212 0.831901    
## as.factor(DID)250       -0.721227   0.105067  -6.864 6.68e-12 ***
## as.factor(DID)251       -0.688719   0.106677  -6.456 1.07e-10 ***
## as.factor(DID)252       -0.913442   0.103706  -8.808  < 2e-16 ***
## as.factor(DID)253       -0.511879   0.111901  -4.574 4.78e-06 ***
## as.factor(DID)254       -1.762981   0.115162 -15.309  < 2e-16 ***
## as.factor(DID)255       -0.856150   0.108863  -7.864 3.71e-15 ***
## as.factor(DID)256       -0.251852   0.106913  -2.356 0.018489 *  
## as.factor(DID)257        0.505590   0.140392   3.601 0.000317 ***
## as.factor(DID)258       -1.562184   0.098968 -15.785  < 2e-16 ***
## as.factor(DID)259       -0.846290   0.107113  -7.901 2.77e-15 ***
## as.factor(DID)260       -0.132057   0.108378  -1.218 0.223037    
## as.factor(DID)261       -2.301424   0.098132 -23.452  < 2e-16 ***
## as.factor(DID)262       -0.725316   0.103528  -7.006 2.45e-12 ***
## as.factor(DID)263       -0.414257   0.107373  -3.858 0.000114 ***
## as.factor(DID)264       -0.444065   0.105615  -4.205 2.62e-05 ***
## as.factor(DID)265       -0.075139   0.114318  -0.657 0.510998    
## as.factor(DID)266        0.615342   0.130073   4.731 2.24e-06 ***
## as.factor(DID)267       -0.257368   0.106773  -2.410 0.015934 *  
## as.factor(DID)268       -1.115655   0.104103 -10.717  < 2e-16 ***
## as.factor(DID)269       -0.143548   0.111944  -1.282 0.199729    
## as.factor(DID)270       -0.052250   0.111599  -0.468 0.639648    
## as.factor(DID)271        0.056815   0.115452   0.492 0.622643    
## as.factor(DID)272       -0.088445   0.112363  -0.787 0.431202    
## as.factor(DID)273       -0.338363   0.109176  -3.099 0.001940 ** 
## as.factor(DID)274       -0.169419   0.120699  -1.404 0.160423    
## as.factor(DID)275       -1.110405   0.101896 -10.897  < 2e-16 ***
## as.factor(DID)276       -0.187979   0.106426  -1.766 0.077348 .  
## as.factor(DID)277       -0.200188   0.146869  -1.363 0.172873    
## as.factor(DID)278       -0.791947   0.098415  -8.047 8.49e-16 ***
## as.factor(DID)279       -1.831521   0.099957 -18.323  < 2e-16 ***
## as.factor(DID)280       -1.045823   0.102210 -10.232  < 2e-16 ***
## as.factor(DID)281       -0.817280   0.101859  -8.024 1.03e-15 ***
## as.factor(DID)282       -0.777911   0.109042  -7.134 9.75e-13 ***
## as.factor(DID)284       -0.782806   0.104522  -7.489 6.92e-14 ***
## as.factor(DID)287       -1.035915   0.104583  -9.905  < 2e-16 ***
## as.factor(DID)289       -0.883532   0.101090  -8.740  < 2e-16 ***
## as.factor(DID)290       -1.049015   0.104241 -10.063  < 2e-16 ***
## as.factor(DID)315       -0.730629   0.115409  -6.331 2.44e-10 ***
## as.factor(DID)316       -1.090704   0.103373 -10.551  < 2e-16 ***
## as.factor(DID)318       -1.666857   0.105254 -15.837  < 2e-16 ***
## as.factor(DID)319       -2.273787   0.098779 -23.019  < 2e-16 ***
## as.factor(DID)320       -1.525834   0.104652 -14.580  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 269656  on 245600  degrees of freedom
## AIC: 269950
## 
## Number of Fisher Scoring iterations: 5
## Warning in ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary =
## TRUE, : maximum three independent variables are allowed
## NULL
Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 415.55, df = 8, p-value < 2.2e-16
exponential transformation
exp(glm_child$coefficients)
##             (Intercept) n_children_in_household          PR004_PR009_01 
##              8.68426215              0.97189458              1.41491756 
##                 C002_01       as.factor(DID)147       as.factor(DID)148 
##              0.51468855              0.54465757              0.91028798 
##       as.factor(DID)149       as.factor(DID)150       as.factor(DID)151 
##              0.44727950              0.45656288              0.82631471 
##       as.factor(DID)152       as.factor(DID)153       as.factor(DID)154 
##              0.71571648              0.40624870              0.50750152 
##       as.factor(DID)155       as.factor(DID)156       as.factor(DID)157 
##              0.43067852              0.26888624              0.63627844 
##       as.factor(DID)158       as.factor(DID)159       as.factor(DID)160 
##              0.70899315              0.26994651              0.47789348 
##       as.factor(DID)161       as.factor(DID)162       as.factor(DID)163 
##              0.27869861              1.15104847              1.07902058 
##       as.factor(DID)164       as.factor(DID)165       as.factor(DID)166 
##              0.56755636              0.37570439              0.43566418 
##       as.factor(DID)167       as.factor(DID)169       as.factor(DID)170 
##              1.66479663              0.17840212              0.62322241 
##       as.factor(DID)171       as.factor(DID)172       as.factor(DID)173 
##              0.64561636              0.52807089              0.76555125 
##       as.factor(DID)174       as.factor(DID)175       as.factor(DID)176 
##              0.38415219              0.26401140              0.62350291 
##       as.factor(DID)177       as.factor(DID)178       as.factor(DID)179 
##              0.42650577              1.22386457              0.48762652 
##       as.factor(DID)180       as.factor(DID)181       as.factor(DID)182 
##              0.50733088              0.88907702              0.24936704 
##       as.factor(DID)183       as.factor(DID)184       as.factor(DID)185 
##              0.51703930              0.27053703              1.33583600 
##       as.factor(DID)186       as.factor(DID)187       as.factor(DID)188 
##              0.18382953              0.23287662              0.27206914 
##       as.factor(DID)189       as.factor(DID)190       as.factor(DID)191 
##              0.69913479              0.37003145              0.40438937 
##       as.factor(DID)192       as.factor(DID)193       as.factor(DID)194 
##              0.34656849              0.94109585              0.18563282 
##       as.factor(DID)195       as.factor(DID)196       as.factor(DID)197 
##              0.15217520              0.49712841              0.23947388 
##       as.factor(DID)198       as.factor(DID)199       as.factor(DID)200 
##              0.16775016              1.42165604              0.23287626 
##       as.factor(DID)202       as.factor(DID)203       as.factor(DID)204 
##              0.22845019              0.23631521              0.40863100 
##       as.factor(DID)205       as.factor(DID)206       as.factor(DID)207 
##              0.20450233              0.19163614              0.13145413 
##       as.factor(DID)208       as.factor(DID)209       as.factor(DID)210 
##              0.22004826              0.35783680              0.17847180 
##       as.factor(DID)211       as.factor(DID)212       as.factor(DID)213 
##              0.35702078              0.26061459              0.20551759 
##       as.factor(DID)214       as.factor(DID)215       as.factor(DID)216 
##              0.12300538              0.25666877              0.18714385 
##       as.factor(DID)217       as.factor(DID)218       as.factor(DID)219 
##              0.18088477              0.19509408              0.21377515 
##       as.factor(DID)220       as.factor(DID)221       as.factor(DID)222 
##              0.12208693              0.08580172              0.18406963 
##       as.factor(DID)223       as.factor(DID)224       as.factor(DID)225 
##              0.20107914              0.11671582              0.15679900 
##       as.factor(DID)226       as.factor(DID)227       as.factor(DID)228 
##              0.24954195              0.18682063              0.10466013 
##       as.factor(DID)229       as.factor(DID)230       as.factor(DID)231 
##              0.13925361              0.11794719              0.24562354 
##       as.factor(DID)232       as.factor(DID)233       as.factor(DID)234 
##              0.74601182              0.10325112              0.20515743 
##       as.factor(DID)235       as.factor(DID)236       as.factor(DID)237 
##              0.43506233              0.33115733              0.46545817 
##       as.factor(DID)238       as.factor(DID)239       as.factor(DID)240 
##              1.28232283              0.40824317              0.28875517 
##       as.factor(DID)241       as.factor(DID)242       as.factor(DID)243 
##              0.09357048              0.33251418              0.34725018 
##       as.factor(DID)244       as.factor(DID)245       as.factor(DID)246 
##              0.57083429              1.29903108              0.18143825 
##       as.factor(DID)247       as.factor(DID)248       as.factor(DID)249 
##              0.29096723              0.38865667              1.02682949 
##       as.factor(DID)250       as.factor(DID)251       as.factor(DID)252 
##              0.48615558              0.50221916              0.40114117 
##       as.factor(DID)253       as.factor(DID)254       as.factor(DID)255 
##              0.59936809              0.17153272              0.42479420 
##       as.factor(DID)256       as.factor(DID)257       as.factor(DID)258 
##              0.77736012              1.65796325              0.20967757 
##       as.factor(DID)259       as.factor(DID)260       as.factor(DID)261 
##              0.42900378              0.87629061              0.10011614 
##       as.factor(DID)262       as.factor(DID)263       as.factor(DID)264 
##              0.48417163              0.66083096              0.64142358 
##       as.factor(DID)265       as.factor(DID)266       as.factor(DID)267 
##              0.92761446              1.85028857              0.77308363 
##       as.factor(DID)268       as.factor(DID)269       as.factor(DID)270 
##              0.32770047              0.86627896              0.94909197 
##       as.factor(DID)271       as.factor(DID)272       as.factor(DID)273 
##              1.05845992              0.91535372              0.71293627 
##       as.factor(DID)274       as.factor(DID)275       as.factor(DID)276 
##              0.84415517              0.32942564              0.82863226 
##       as.factor(DID)277       as.factor(DID)278       as.factor(DID)279 
##              0.81857708              0.45296190              0.16016976 
##       as.factor(DID)280       as.factor(DID)281       as.factor(DID)282 
##              0.35140253              0.44163131              0.45936450 
##       as.factor(DID)284       as.factor(DID)287       as.factor(DID)289 
##              0.45712131              0.35490148              0.41332046 
##       as.factor(DID)290       as.factor(DID)315       as.factor(DID)316 
##              0.35028265              0.48160583              0.33598000 
##       as.factor(DID)318       as.factor(DID)319       as.factor(DID)320 
##              0.18883967              0.10292172              0.21743960
confidence interval (intercept and coefficient)
# confint(glm_child, level = 0.95)
exponential transformation of confidence interval (odds ratio of intercept and coefficient)
# exp(confint(glm_child, level = 0.95))
AIC
extractAIC(glm_child)
## [1]    147.0 269950.4
BIC
extractAIC(glm_child, k = log(nrow(glm_child$data)))
## [1]    147 271481
effectiveness of explanatory variables
glm_child_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy)
anova(glm_child_null, glm_child, test = "Chisq")
variables selection
# step(glm_child_null, direction = "both", 
#      scope = (~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(DID)))
multicolinearity
vif(glm_child)
##                             GVIF  Df GVIF^(1/(2*Df))
## n_children_in_household 1.146849   1        1.070910
## PR004_PR009_01          1.246852   1        1.116625
## C002_01                 1.024292   1        1.012073
## as.factor(DID)          1.441880 143        1.001280
Districts with P-Values Bigger Than 0.05 and Coefficient
glm_summary <- glm_child %>% summary()
glm_coef <- as.data.frame(glm_summary$coefficients)
glm_coef %>% filter(`Pr(>|z|)` > 0.05)
Districts with P-Values Bigger Than 0.05
glm_coef %>% filter(`Pr(>|z|)` > 0.05) %>% t() %>% .[-c(1:5),]
##      as.factor(DID)148 as.factor(DID)151 as.factor(DID)162 as.factor(DID)163
##      as.factor(DID)178 as.factor(DID)181 as.factor(DID)193 as.factor(DID)249
##      as.factor(DID)260 as.factor(DID)265 as.factor(DID)269 as.factor(DID)270
##      as.factor(DID)271 as.factor(DID)272 as.factor(DID)274 as.factor(DID)276
##      as.factor(DID)277
dists_not_fit <-  c(148, 151, 162, 163, 178, 181, 193, 249, 260, 265, 269, 270, 271, 272, 274, 276, 277)
child_ica_dummy %>% filter(DID == dists_not_fit) %>% select(DID, DNAME) %>% summarize(DID = unique(DID), DNAME = unique(DNAME))
## Warning in DID == dists_not_fit: 長いオブジェクトの長さが短いオブジェクトの長さ
## の倍数になっていません
Districts with P-Values Bigger Than 0.001 and Coefficient
glm_summary <- glm_child %>% summary()
glm_coef <- as.data.frame(glm_summary$coefficients)
glm_coef %>% filter(`Pr(>|z|)` > 0.001)
Districts with P-Values Bigger Than 0.001
glm_coef %>% filter(`Pr(>|z|)` > 0.001) %>% t() %>% .[-c(1:5),]
##      as.factor(DID)148 as.factor(DID)151 as.factor(DID)152 as.factor(DID)158
##      as.factor(DID)162 as.factor(DID)163 as.factor(DID)173 as.factor(DID)178
##      as.factor(DID)181 as.factor(DID)185 as.factor(DID)193 as.factor(DID)199
##      as.factor(DID)232 as.factor(DID)238 as.factor(DID)245 as.factor(DID)249
##      as.factor(DID)256 as.factor(DID)260 as.factor(DID)265 as.factor(DID)267
##      as.factor(DID)269 as.factor(DID)270 as.factor(DID)271 as.factor(DID)272
##      as.factor(DID)273 as.factor(DID)274 as.factor(DID)276 as.factor(DID)277
dists_not_fit <-  c(148, 151, 152, 158, 162, 163, 173, 178, 181, 185, 193, 199, 232, 238, 245, 249, 256, 260, 265, 267, 269, 270, 271, 272, 273, 274, 276, 277)
child_ica_dummy %>% filter(DID == dists_not_fit) %>% select(DID, DNAME) %>% summarize(DID = unique(DID), DNAME = unique(DNAME))
## Warning in DID == dists_not_fit: 長いオブジェクトの長さが短いオブジェクトの長さ
## の倍数になっていません

Dists_not_fit. GLM C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01

All of dists_not_fit
glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(VID), family = "binomial", data = child_ica_dummy %>% filter(DID %in% dists_not_fit))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01 + as.factor(VID), family = "binomial", data = child_ica_dummy %>% 
##     filter(DID %in% dists_not_fit))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.9975   0.2123   0.4372   0.5768   1.4738  
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               1.702857   0.485396   3.508 0.000451 ***
## n_children_in_household   0.036001   0.012295   2.928 0.003410 ** 
## PR004_PR009_01            0.085176   0.040443   2.106 0.035197 *  
## C002_01                  -0.181653   0.028957  -6.273 3.54e-10 ***
## as.factor(VID)4937       -0.475516   0.669974  -0.710 0.477857    
## as.factor(VID)4938        1.793771   1.123676   1.596 0.110412    
## as.factor(VID)4939        1.106393   0.872792   1.268 0.204924    
## as.factor(VID)4940        0.394658   0.715379   0.552 0.581170    
## as.factor(VID)4941       -0.345331   0.662635  -0.521 0.602264    
## as.factor(VID)4942        1.897602   1.122701   1.690 0.090987 .  
## as.factor(VID)4943       -1.210639   0.566444  -2.137 0.032577 *  
## as.factor(VID)4944        0.575980   0.774371   0.744 0.456995    
## as.factor(VID)4945        0.078486   0.614769   0.128 0.898411    
## as.factor(VID)4946        0.527867   0.712682   0.741 0.458890    
## as.factor(VID)4947       -0.398360   0.623935  -0.638 0.523172    
## as.factor(VID)4948        1.326896   0.869342   1.526 0.126929    
## as.factor(VID)4949       -0.512888   0.613484  -0.836 0.403141    
## as.factor(VID)4950        2.151462   1.119483   1.922 0.054627 .  
## as.factor(VID)4951        0.858204   0.768481   1.117 0.264100    
## as.factor(VID)4952       -0.381300   0.580474  -0.657 0.511260    
## as.factor(VID)4953        1.603477   0.866352   1.851 0.064193 .  
## as.factor(VID)4954       -1.093095   0.577997  -1.891 0.058601 .  
## as.factor(VID)4955       -0.050238   0.633927  -0.079 0.936835    
## as.factor(VID)4956       -1.544766   0.538628  -2.868 0.004131 ** 
## as.factor(VID)4957        0.048790   0.682060   0.072 0.942973    
## as.factor(VID)4958       -0.811464   0.565449  -1.435 0.151264    
## as.factor(VID)4959       15.776122 688.099301   0.023 0.981708    
## as.factor(VID)4966       -0.047350   0.616762  -0.077 0.938805    
## as.factor(VID)4991       15.732527 482.918751   0.033 0.974011    
## as.factor(VID)4992        0.927659   0.767274   1.209 0.226650    
## as.factor(VID)4993        0.508508   0.775083   0.656 0.511780    
## as.factor(VID)4994        2.147531   1.119164   1.919 0.055001 .  
## as.factor(VID)4995        0.475163   0.713868   0.666 0.505656    
## as.factor(VID)4996        0.240871   0.678371   0.355 0.722535    
## as.factor(VID)4997        1.789771   1.124045   1.592 0.111326    
## as.factor(VID)4998        0.468174   0.775809   0.603 0.546200    
## as.factor(VID)4999        1.570511   0.866595   1.812 0.069943 .  
## as.factor(VID)5000        0.492290   0.712538   0.691 0.489631    
## as.factor(VID)5001        0.661238   0.771970   0.857 0.391689    
## as.factor(VID)5002       -0.341784   0.641394  -0.533 0.594119    
## as.factor(VID)14940       0.561585   0.712159   0.789 0.430365    
## as.factor(VID)14941       0.316169   0.647988   0.488 0.625604    
## as.factor(VID)14942      -0.630604   0.604462  -1.043 0.296833    
## as.factor(VID)14943       0.467471   0.674064   0.694 0.487989    
## as.factor(VID)14944       0.461904   0.674483   0.685 0.493453    
## as.factor(VID)14945       2.049224   1.120379   1.829 0.067393 .  
## as.factor(VID)14947       2.057430   1.120162   1.837 0.066250 .  
## as.factor(VID)14950       0.296842   0.648604   0.458 0.647194    
## as.factor(VID)14953      -0.062520   0.634450  -0.099 0.921501    
## as.factor(VID)14954      -0.374583   0.588340  -0.637 0.524335    
## as.factor(VID)14955       0.327471   0.676200   0.484 0.628186    
## as.factor(VID)14960       1.302113   0.869561   1.497 0.134279    
## as.factor(VID)14961       1.189938   0.871560   1.365 0.172160    
## as.factor(VID)14963       0.255933   0.677879   0.378 0.705765    
## as.factor(VID)14968       0.908239   0.767804   1.183 0.236847    
## as.factor(VID)14974      15.765453 602.859981   0.026 0.979137    
## as.factor(VID)14989      -1.924316   0.557627  -3.451 0.000559 ***
## as.factor(VID)14992       0.626077   0.712945   0.878 0.379859    
## as.factor(VID)15044       0.031497   0.615695   0.051 0.959201    
## as.factor(VID)15045       0.507392   0.713008   0.712 0.476699    
## as.factor(VID)15046      -0.193405   0.605654  -0.319 0.749475    
## as.factor(VID)15047       1.124288   0.764798   1.470 0.141549    
## as.factor(VID)15048      15.738471 475.904668   0.033 0.973618    
## as.factor(VID)15049       0.239993   0.649539   0.369 0.711768    
## as.factor(VID)15050       1.059699   0.765894   1.384 0.166478    
## as.factor(VID)15051       0.087376   0.614193   0.142 0.886873    
## as.factor(VID)15052       2.162975   1.119121   1.933 0.053268 .  
## as.factor(VID)15053       1.549129   0.866415   1.788 0.073780 .  
## as.factor(VID)15054       2.442462   1.117262   2.186 0.028807 *  
## as.factor(VID)15055      15.760623 465.867725   0.034 0.973012    
## as.factor(VID)15056       0.016514   0.653380   0.025 0.979836    
## as.factor(VID)15057      -0.311420   0.608224  -0.512 0.608641    
## as.factor(VID)15058       0.560069   0.672762   0.832 0.405131    
## as.factor(VID)15059      -0.745893   0.581300  -1.283 0.199441    
## as.factor(VID)15060      -0.317363   0.598154  -0.531 0.595716    
## as.factor(VID)15061       0.087698   0.601177   0.146 0.884018    
## as.factor(VID)15062       0.168002   0.629668   0.267 0.789615    
## as.factor(VID)15063      -0.595903   0.577550  -1.032 0.302176    
## as.factor(VID)15064       0.343953   0.626848   0.549 0.583210    
## as.factor(VID)15065       1.323116   0.762586   1.735 0.082734 .  
## as.factor(VID)15066      -0.058239   0.603119  -0.097 0.923073    
## as.factor(VID)15067      -0.085077   0.604256  -0.141 0.888032    
## as.factor(VID)15068       0.447159   0.674373   0.663 0.507284    
## as.factor(VID)15069      -0.621296   0.593878  -1.046 0.295484    
## as.factor(VID)15071      -0.217174   0.595648  -0.365 0.715409    
## as.factor(VID)15072      -0.707899   0.588761  -1.202 0.229227    
## as.factor(VID)15073       0.336505   0.648069   0.519 0.603592    
## as.factor(VID)15074      -0.156413   0.593777  -0.263 0.792226    
## as.factor(VID)15075       2.215300   1.118873   1.980 0.047710 *  
## as.factor(VID)15076      -0.500689   0.582727  -0.859 0.390220    
## as.factor(VID)15077      -0.725065   0.574332  -1.262 0.206787    
## as.factor(VID)15078       0.031167   0.614725   0.051 0.959565    
## as.factor(VID)15079      -0.994810   0.564925  -1.761 0.078245 .  
## as.factor(VID)15080       0.780988   0.708754   1.102 0.270498    
## as.factor(VID)15081      -0.492756   0.569364  -0.865 0.386792    
## as.factor(VID)15082       0.035333   0.615607   0.057 0.954230    
## as.factor(VID)15083       0.904192   0.707535   1.278 0.201268    
## as.factor(VID)15084      -0.178550   0.584903  -0.305 0.760165    
## as.factor(VID)15085      -0.076212   0.603747  -0.126 0.899548    
## as.factor(VID)15086      -0.148628   0.605298  -0.246 0.806034    
## as.factor(VID)15087       2.279669   1.118420   2.038 0.041521 *  
## as.factor(VID)15088       0.525381   0.645341   0.814 0.415580    
## as.factor(VID)15089       0.202400   0.629367   0.322 0.747760    
## as.factor(VID)15090       0.594010   0.644599   0.922 0.356780    
## as.factor(VID)15091      -0.153131   0.593783  -0.258 0.796492    
## as.factor(VID)15092       0.100610   0.614164   0.164 0.869876    
## as.factor(VID)15093       0.396640   0.626466   0.633 0.526644    
## as.factor(VID)15094       2.209155   1.118887   1.974 0.048334 *  
## as.factor(VID)15095       0.454650   0.625564   0.727 0.467358    
## as.factor(VID)15096       0.544041   0.645180   0.843 0.399095    
## as.factor(VID)15097       0.798799   0.669538   1.193 0.232845    
## as.factor(VID)15098      15.713542 459.537358   0.034 0.972722    
## as.factor(VID)15099      -0.257762   0.556321  -0.463 0.643125    
## as.factor(VID)15100      -0.111448   0.634432  -0.176 0.860556    
## as.factor(VID)15101       0.828043   0.708398   1.169 0.242446    
## as.factor(VID)15102      -0.759386   0.552564  -1.374 0.169350    
## as.factor(VID)15103       0.933183   0.708297   1.318 0.187671    
## as.factor(VID)15104      15.774111 602.845083   0.026 0.979125    
## as.factor(VID)15105       1.158626   0.764505   1.516 0.129640    
## as.factor(VID)15106      -1.541444   0.581038  -2.653 0.007980 ** 
## as.factor(VID)15107       0.092825   0.601300   0.154 0.877315    
## as.factor(VID)15108       0.783610   0.671117   1.168 0.242960    
## as.factor(VID)15109      -0.275224   0.579003  -0.475 0.634544    
## as.factor(VID)15110      -1.097041   0.552123  -1.987 0.046928 *  
## as.factor(VID)15111       0.421882   0.647115   0.652 0.514438    
## as.factor(VID)15112       1.472668   0.867615   1.697 0.089626 .  
## as.factor(VID)15113       0.601049   0.672483   0.894 0.371442    
## as.factor(VID)15114       0.058866   0.601160   0.098 0.921995    
## as.factor(VID)15115      -0.949906   0.546634  -1.738 0.082257 .  
## as.factor(VID)15116      -0.884227   0.585786  -1.509 0.131179    
## as.factor(VID)15117       0.185356   0.680132   0.273 0.785215    
## as.factor(VID)15118       1.176895   0.764287   1.540 0.123594    
## as.factor(VID)15119       1.317172   0.764104   1.724 0.084742 .  
## as.factor(VID)15120       0.979764   0.707162   1.385 0.165904    
## as.factor(VID)15121      -0.098839   0.605019  -0.163 0.870231    
## as.factor(VID)15122       0.187665   0.629900   0.298 0.765757    
## as.factor(VID)15123       1.697110   0.865835   1.960 0.049986 *  
## as.factor(VID)15124       2.456604   1.117341   2.199 0.027905 *  
## as.factor(VID)15125       0.575170   0.672692   0.855 0.392536    
## as.factor(VID)15126      -0.937062   0.546165  -1.716 0.086215 .  
## as.factor(VID)15127       0.249629   0.628502   0.397 0.691234    
## as.factor(VID)15128       0.141376   0.601193   0.235 0.814085    
## as.factor(VID)15129       0.090168   0.590226   0.153 0.878580    
## as.factor(VID)15130       0.379718   0.648300   0.586 0.558068    
## as.factor(VID)15131       0.568973   0.672977   0.845 0.397856    
## as.factor(VID)15132      -1.145392   0.574300  -1.994 0.046107 *  
## as.factor(VID)15133      -1.269375   0.597001  -2.126 0.033482 *  
## as.factor(VID)15134      -0.088010   0.633871  -0.139 0.889572    
## as.factor(VID)15135      -0.713519   0.588092  -1.213 0.225023    
## as.factor(VID)15136      -1.166130   0.556618  -2.095 0.036168 *  
## as.factor(VID)15137      -0.803828   0.581920  -1.381 0.167175    
## as.factor(VID)15138      -0.970739   0.569465  -1.705 0.088260 .  
## as.factor(VID)15139       0.594900   0.672559   0.885 0.376409    
## as.factor(VID)15140       0.379461   0.647442   0.586 0.557814    
## as.factor(VID)15141       0.211660   0.628971   0.337 0.736480    
## as.factor(VID)15142       2.229658   1.118739   1.993 0.046260 *  
## as.factor(VID)15143      -0.817939   0.564555  -1.449 0.147387    
## as.factor(VID)15144       0.496236   0.645615   0.769 0.442116    
## as.factor(VID)15145       1.060963   0.705748   1.503 0.132757    
## as.factor(VID)15168      -0.875645   0.554252  -1.580 0.114137    
## as.factor(VID)15171       0.444003   0.776587   0.572 0.567501    
## as.factor(VID)15173      -0.397876   0.567723  -0.701 0.483411    
## as.factor(VID)15179      -0.784681   0.569580  -1.378 0.168312    
## as.factor(VID)15187      -1.108795   0.599173  -1.851 0.064235 .  
## as.factor(VID)15190       0.195343   0.651127   0.300 0.764171    
## as.factor(VID)15192      -0.234835   0.595301  -0.394 0.693226    
## as.factor(VID)15195       0.277659   0.677956   0.410 0.682133    
## as.factor(VID)15197       1.179884   0.871594   1.354 0.175829    
## as.factor(VID)15198      15.812142 747.168247   0.021 0.983116    
## as.factor(VID)15201       1.125265   0.872770   1.289 0.197293    
## as.factor(VID)15202      -0.443955   0.581914  -0.763 0.445510    
## as.factor(VID)15205       0.720569   0.770940   0.935 0.349962    
## as.factor(VID)15207       0.214382   0.650666   0.329 0.741792    
## as.factor(VID)15210      -0.607759   0.565675  -1.074 0.282645    
## as.factor(VID)15211       0.713257   0.710198   1.004 0.315231    
## as.factor(VID)15213       0.120452   0.651953   0.185 0.853420    
## as.factor(VID)15229       0.776140   0.709195   1.094 0.273782    
## as.factor(VID)15230      -0.184525   0.606469  -0.304 0.760929    
## as.factor(VID)15232      15.765035 559.073195   0.028 0.977504    
## as.factor(VID)15235      -0.021270   0.632560  -0.034 0.973175    
## as.factor(VID)15236       0.414324   0.675056   0.614 0.539372    
## as.factor(VID)15239      -0.391143   0.641310  -0.610 0.541919    
## as.factor(VID)15241       1.055754   0.765591   1.379 0.167893    
## as.factor(VID)15243      15.783738 632.976384   0.025 0.980106    
## as.factor(VID)15244       0.442885   0.714340   0.620 0.535263    
## as.factor(VID)15247      -0.448461   0.589688  -0.761 0.446952    
## as.factor(VID)15249      -0.801271   0.620566  -1.291 0.196636    
## as.factor(VID)15251       0.812615   0.769214   1.056 0.290775    
## as.factor(VID)15254      -0.919163   0.578628  -1.589 0.112168    
## as.factor(VID)15255       0.109642   0.631475   0.174 0.862158    
## as.factor(VID)15257      -0.450750   0.626377  -0.720 0.471762    
## as.factor(VID)15260      -0.484237   0.600305  -0.807 0.419867    
## as.factor(VID)15261       2.019022   1.120569   1.802 0.071579 .  
## as.factor(VID)15262       0.312191   0.717129   0.435 0.663320    
## as.factor(VID)15263       1.224011   0.871578   1.404 0.160211    
## as.factor(VID)15264       0.354613   0.676181   0.524 0.599976    
## as.factor(VID)15265      -0.738401   0.588066  -1.256 0.209245    
## as.factor(VID)15267      -0.784789   0.560254  -1.401 0.161281    
## as.factor(VID)15268      -1.032776   0.606968  -1.702 0.088843 .  
## as.factor(VID)15269      -0.545702   0.584067  -0.934 0.350142    
## as.factor(VID)15270      -1.360080   0.600804  -2.264 0.023588 *  
## as.factor(VID)15272      -0.481800   0.699814  -0.688 0.491158    
## as.factor(VID)15276      -0.004237   0.654524  -0.006 0.994835    
## as.factor(VID)15278      -0.560398   0.629330  -0.890 0.373215    
## as.factor(VID)15279      -0.014631   0.632964  -0.023 0.981558    
## as.factor(VID)15280      -0.546942   0.602742  -0.907 0.364184    
## as.factor(VID)15285      -0.215053   0.620202  -0.347 0.728781    
## as.factor(VID)15286       0.478304   0.714303   0.670 0.503107    
## as.factor(VID)15288      15.797358 633.142810   0.025 0.980094    
## as.factor(VID)15289      -0.013204   0.602969  -0.022 0.982529    
## as.factor(VID)15293      15.771705 595.972267   0.026 0.978887    
## as.factor(VID)15349      -0.894303   0.613082  -1.459 0.144648    
## as.factor(VID)15360      -0.678722   0.618161  -1.098 0.272218    
## as.factor(VID)15403       1.979355   1.120986   1.766 0.077442 .  
## as.factor(VID)15404      15.732195 479.415460   0.033 0.973822    
## as.factor(VID)15407       0.536933   0.672961   0.798 0.424948    
## as.factor(VID)15409       0.710553   0.670653   1.059 0.289374    
## as.factor(VID)15413      -0.395649   0.598063  -0.662 0.508259    
## as.factor(VID)15414       0.137504   0.613914   0.224 0.822773    
## as.factor(VID)15422       0.367269   0.676198   0.543 0.587035    
## as.factor(VID)15425      -0.069733   0.655686  -0.106 0.915303    
## as.factor(VID)15427       0.344312   0.647818   0.531 0.595076    
## as.factor(VID)15432      15.762841 528.337912   0.030 0.976199    
## as.factor(VID)15435      -0.673387   0.558023  -1.207 0.227534    
## as.factor(VID)15436       0.045296   0.632108   0.072 0.942873    
## as.factor(VID)15437      -0.838054   0.590380  -1.420 0.155749    
## as.factor(VID)15445      -1.854087   0.549442  -3.374 0.000740 ***
## as.factor(VID)25434       0.585411   0.712241   0.822 0.411118    
## as.factor(VID)25435      -0.294754   0.596412  -0.494 0.621157    
## as.factor(VID)25441      -0.606496   0.593873  -1.021 0.307134    
## as.factor(VID)25442       0.645421   0.644100   1.002 0.316319    
## as.factor(VID)25443       0.124661   0.630942   0.198 0.843375    
## as.factor(VID)25447       0.241576   0.628328   0.384 0.700627    
## as.factor(VID)25450       0.737328   0.709115   1.040 0.298439    
## as.factor(VID)25452       0.045977   0.653371   0.070 0.943900    
## as.factor(VID)25455      -1.196015   0.561468  -2.130 0.033159 *  
## as.factor(VID)25456      -0.324135   0.561735  -0.577 0.563922    
## as.factor(VID)35435      -0.450432   0.568517  -0.792 0.428190    
## as.factor(VID)35439      15.727720 514.677160   0.031 0.975622    
## as.factor(VID)35440       0.617370   0.710828   0.869 0.385108    
## as.factor(VID)35447       0.518947   0.645625   0.804 0.421518    
## as.factor(VID)35448       0.604583   0.672021   0.900 0.368307    
## as.factor(VID)35449       0.868802   0.707955   1.227 0.219747    
## as.factor(VID)35456       1.070299   0.765691   1.398 0.162167    
## as.factor(VID)35457      -0.801998   0.569609  -1.408 0.159137    
## as.factor(VID)35460       0.731820   0.709349   1.032 0.302223    
## as.factor(VID)35462      -0.808360   0.557586  -1.450 0.147128    
## as.factor(VID)35467       0.208587   0.612675   0.340 0.733516    
## as.factor(VID)45434       0.232671   0.599885   0.388 0.698121    
## as.factor(VID)45436       0.525345   0.673183   0.780 0.435161    
## as.factor(VID)45441       1.524237   0.866902   1.758 0.078704 .  
## as.factor(VID)45443       2.324441   1.117762   2.080 0.037567 *  
## as.factor(VID)45444      -0.524629   0.555738  -0.944 0.345158    
## as.factor(VID)45447       0.536349   0.645410   0.831 0.405962    
## as.factor(VID)45448       0.268813   0.628307   0.428 0.668770    
## as.factor(VID)45451      -0.472789   0.599813  -0.788 0.430564    
## as.factor(VID)45452       0.130439   0.613682   0.213 0.831677    
## as.factor(VID)45453      -0.137969   0.594074  -0.232 0.816350    
## as.factor(VID)45454      -0.047317   0.604363  -0.078 0.937596    
## as.factor(VID)45455      15.716890 441.875937   0.036 0.971626    
## as.factor(VID)45456      -0.474591   0.582718  -0.814 0.415391    
## as.factor(VID)45457      -1.546024   0.561348  -2.754 0.005885 ** 
## as.factor(VID)45458      -0.229574   0.595667  -0.385 0.699936    
## as.factor(VID)45461      -0.656476   0.586178  -1.120 0.262746    
## as.factor(VID)45462      -1.217592   0.566330  -2.150 0.031558 *  
## as.factor(VID)45463       1.283525   0.763573   1.681 0.092773 .  
## as.factor(VID)45464      -0.196526   0.606629  -0.324 0.745965    
## as.factor(VID)45469       2.258281   1.118593   2.019 0.043502 *  
## as.factor(VID)45470       0.245368   0.649480   0.378 0.705586    
## as.factor(VID)45472       1.382932   0.868538   1.592 0.111328    
## as.factor(VID)45473       0.181447   0.600406   0.302 0.762494    
## as.factor(VID)45474       1.284656   0.870479   1.476 0.139996    
## as.factor(VID)45476      -1.071723   0.571650  -1.875 0.060822 .  
## as.factor(VID)45480      -0.743074   0.574249  -1.294 0.195668    
## as.factor(VID)45481       0.283633   0.677036   0.419 0.675265    
## as.factor(VID)45483       1.072220   0.765722   1.400 0.161432    
## as.factor(VID)45485       0.871936   0.708430   1.231 0.218398    
## as.factor(VID)45486       0.172235   0.651434   0.264 0.791477    
## as.factor(VID)45487      -0.445029   0.612739  -0.726 0.467658    
## as.factor(VID)45500      -0.138878   0.636102  -0.218 0.827175    
## as.factor(VID)45501       0.726732   0.709875   1.024 0.305955    
## as.factor(VID)45502      -0.056545   0.633424  -0.089 0.928869    
## as.factor(VID)45504      -0.417609   0.611194  -0.683 0.494438    
## as.factor(VID)45506       0.677067   0.710472   0.953 0.340599    
## as.factor(VID)45508       0.520706   0.673713   0.773 0.439587    
## as.factor(VID)45511       0.165623   0.613028   0.270 0.787028    
## as.factor(VID)45513       0.529866   0.645132   0.821 0.411458    
## as.factor(VID)45514      -0.299359   0.572856  -0.523 0.601271    
## as.factor(VID)45517       0.831440   0.769178   1.081 0.279721    
## as.factor(VID)45519      15.774323 548.189197   0.029 0.977044    
## as.factor(VID)45522      15.732402 510.246410   0.031 0.975403    
## as.factor(VID)45523       1.493238   0.867513   1.721 0.085199 .  
## as.factor(VID)45526      -1.132489   0.557153  -2.033 0.042089 *  
## as.factor(VID)45527       0.100948   0.652817   0.155 0.877109    
## as.factor(VID)45528      -0.425574   0.563454  -0.755 0.450072    
## as.factor(VID)45532       0.185870   0.629179   0.295 0.767676    
## as.factor(VID)45533      -0.268175   0.579253  -0.463 0.643388    
## as.factor(VID)45536      -0.388042   0.580937  -0.668 0.504160    
## as.factor(VID)45539      15.750972 510.248825   0.031 0.975374    
## as.factor(VID)45541      -0.782507   0.556060  -1.407 0.159358    
## as.factor(VID)45545      -0.332823   0.597184  -0.557 0.577308    
## as.factor(VID)45547      15.761490 641.378408   0.025 0.980394    
## as.factor(VID)45552      15.743618 542.830047   0.029 0.976862    
## as.factor(VID)45553      -0.015536   0.616364  -0.025 0.979891    
## as.factor(VID)45554      -0.215555   0.585830  -0.368 0.712911    
## as.factor(VID)45555       0.495231   0.646543   0.766 0.443696    
## as.factor(VID)45557      -0.410941   0.599169  -0.686 0.492807    
## as.factor(VID)45558       0.380863   0.626879   0.608 0.543483    
## as.factor(VID)45559      -2.330282   0.601153  -3.876 0.000106 ***
## as.factor(VID)45560       0.061174   0.614594   0.100 0.920713    
## as.factor(VID)45564      -0.822045   0.561263  -1.465 0.143021    
## as.factor(VID)45568       0.026803   0.615699   0.044 0.965276    
## as.factor(VID)45569      -0.514605   0.612627  -0.840 0.400910    
## as.factor(VID)45571      -0.137921   0.593837  -0.232 0.816341    
## as.factor(VID)45574       0.979173   0.874842   1.119 0.263031    
## as.factor(VID)45577      -0.140681   0.594044  -0.237 0.812797    
## as.factor(VID)45579      -1.634751   0.548132  -2.982 0.002860 ** 
## as.factor(VID)45581      -0.825719   0.556801  -1.483 0.138083    
## as.factor(VID)45582       0.513828   0.774981   0.663 0.507318    
## as.factor(VID)45586      -0.881751   0.555400  -1.588 0.112378    
## as.factor(VID)45587       1.129966   0.764931   1.477 0.139619    
## as.factor(VID)45594       0.587489   0.672435   0.874 0.382296    
## as.factor(VID)45596       0.220667   0.600425   0.368 0.713232    
## as.factor(VID)45599      -0.252534   0.691137  -0.365 0.714821    
## as.factor(VID)45601      -1.059780   0.550002  -1.927 0.053996 .  
## as.factor(VID)45602       0.456311   0.714312   0.639 0.522945    
## as.factor(VID)45605      15.835380 734.172271   0.022 0.982792    
## as.factor(VID)45609      15.717901 519.108460   0.030 0.975845    
## as.factor(VID)45610       0.166376   0.614403   0.271 0.786551    
## as.factor(VID)45612      -0.193026   0.619634  -0.312 0.755408    
## as.factor(VID)45614       0.027462   0.654525   0.042 0.966532    
## as.factor(VID)45615       0.162436   0.720961   0.225 0.821742    
## as.factor(VID)45617      15.712402 453.376895   0.035 0.972354    
## as.factor(VID)45622       1.405669   0.868765   1.618 0.105661    
## as.factor(VID)45624      15.735830 465.733436   0.034 0.973047    
## as.factor(VID)45625      -0.248969   0.595795  -0.418 0.676037    
## as.factor(VID)45627       0.967725   0.874755   1.106 0.268605    
## as.factor(VID)45629      -0.270400   0.566900  -0.477 0.633376    
## as.factor(VID)45634       0.854179   0.707684   1.207 0.227430    
## as.factor(VID)45641      -0.655566   0.561982  -1.167 0.243402    
## as.factor(VID)45642       0.571762   0.672797   0.850 0.395420    
## as.factor(VID)45644      -0.096883   0.570258  -0.170 0.865093    
## as.factor(VID)45651       0.041471   0.591515   0.070 0.944106    
## as.factor(VID)45652      15.789574 570.607939   0.028 0.977924    
## as.factor(VID)45653      -0.746596   0.563669  -1.325 0.185327    
## as.factor(VID)45655       0.364385   0.627003   0.581 0.561136    
## as.factor(VID)45661       0.561570   0.673187   0.834 0.404170    
## as.factor(VID)45667      -0.497560   0.570056  -0.873 0.382758    
## as.factor(VID)45672       0.560913   0.645402   0.869 0.384798    
## as.factor(VID)45677       0.304330   0.648554   0.469 0.638895    
## as.factor(VID)45678       0.007467   0.602734   0.012 0.990115    
## as.factor(VID)45683       0.983266   0.766862   1.282 0.199775    
## as.factor(VID)45686       0.774620   0.770075   1.006 0.314462    
## as.factor(VID)45694      -0.093487   0.617813  -0.151 0.879724    
## as.factor(VID)45698      -0.242125   0.586359  -0.413 0.679658    
## as.factor(VID)45702       2.149883   1.119482   1.920 0.054804 .  
## as.factor(VID)45705      -0.161485   0.605154  -0.267 0.789585    
## as.factor(VID)45707       1.075386   0.765450   1.405 0.160049    
## as.factor(VID)45711      -0.335596   0.608850  -0.551 0.581498    
## as.factor(VID)45723       2.060077   1.120396   1.839 0.065959 .  
## as.factor(VID)45734      -0.148588   0.594400  -0.250 0.802603    
## as.factor(VID)45735      -1.498172   0.553934  -2.705 0.006839 ** 
## as.factor(VID)45738       0.135547   0.601065   0.226 0.821581    
## as.factor(VID)45741      -0.171082   0.585728  -0.292 0.770223    
## as.factor(VID)45748      -0.778700   0.564463  -1.380 0.167728    
## as.factor(VID)45758       1.736849   0.865089   2.008 0.044674 *  
## as.factor(VID)45770       0.084249   0.601471   0.140 0.888604    
## as.factor(VID)45772       0.439703   0.675096   0.651 0.514840    
## as.factor(VID)45773      -0.124482   0.604915  -0.206 0.836959    
## as.factor(VID)45775      -0.142935   0.634369  -0.225 0.821732    
## as.factor(VID)45776       0.329473   0.716623   0.460 0.645690    
## as.factor(VID)45778       0.357285   0.626768   0.570 0.568648    
## as.factor(VID)45779      -0.892544   0.584110  -1.528 0.126502    
## as.factor(VID)45781      -0.533461   0.583562  -0.914 0.360640    
## as.factor(VID)45784       0.267886   0.612136   0.438 0.661658    
## as.factor(VID)45785       0.286529   0.718449   0.399 0.690028    
## as.factor(VID)45786       0.844231   0.878175   0.961 0.336378    
## as.factor(VID)45787       0.153146   0.600379   0.255 0.798659    
## as.factor(VID)45788       1.701730   0.865618   1.966 0.049309 *  
## as.factor(VID)45789      15.744107 519.034778   0.030 0.975801    
## as.factor(VID)45791       0.242094   0.612591   0.395 0.692697    
## as.factor(VID)45792      -0.351566   0.597573  -0.588 0.556316    
## as.factor(VID)45793      -0.152159   0.593986  -0.256 0.797823    
## as.factor(VID)45794      -0.850461   0.570679  -1.490 0.136156    
## as.factor(VID)45795      -0.555695   0.601490  -0.924 0.355558    
## as.factor(VID)45796       0.097605   0.614618   0.159 0.873822    
## as.factor(VID)45798      -0.330491   0.662659  -0.499 0.617966    
## as.factor(VID)45800       0.853856   0.708215   1.206 0.227954    
## as.factor(VID)45801      -0.077748   0.603515  -0.129 0.897496    
## as.factor(VID)45802       1.238645   0.870662   1.423 0.154838    
## as.factor(VID)45804       1.966173   1.120382   1.755 0.079274 .  
## as.factor(VID)45806       0.081086   0.615237   0.132 0.895145    
## as.factor(VID)45810      -0.209514   0.578545  -0.362 0.717248    
## as.factor(VID)45815      15.708415 456.352187   0.034 0.972541    
## as.factor(VID)45818      -0.153379   0.593984  -0.258 0.796236    
## as.factor(VID)45823       1.561843   0.866594   1.802 0.071502 .  
## as.factor(VID)45826       1.254044   0.870129   1.441 0.149524    
## as.factor(VID)45828       0.724834   0.709668   1.021 0.307079    
## as.factor(VID)45836      -0.448371   0.574978  -0.780 0.435506    
## as.factor(VID)45838      -0.362397   0.580317  -0.624 0.532312    
## as.factor(VID)45845      -0.259078   0.586969  -0.441 0.658935    
## as.factor(VID)45848      -0.317857   0.587606  -0.541 0.588552    
## as.factor(VID)45854      -0.240073   0.578850  -0.415 0.678332    
## as.factor(VID)45868       0.701155   0.670718   1.045 0.295848    
## as.factor(VID)45894      -1.111706   0.542136  -2.051 0.040306 *  
## as.factor(VID)45907       0.382541   0.626584   0.611 0.541518    
## as.factor(VID)45915       0.629605   0.644144   0.977 0.328357    
## as.factor(VID)45953       1.470617   0.868264   1.694 0.090314 .  
## as.factor(VID)45963      -0.659903   0.595371  -1.108 0.267694    
## as.factor(VID)45967       1.075822   0.765383   1.406 0.159843    
## as.factor(VID)45968      -1.084220   0.576236  -1.882 0.059896 .  
## as.factor(VID)45972       1.621541   0.865892   1.873 0.061112 .  
## as.factor(VID)45978      -0.128046   0.576888  -0.222 0.824345    
## as.factor(VID)45981       1.675045   0.865350   1.936 0.052906 .  
## as.factor(VID)45982      -0.457264   0.582316  -0.785 0.432306    
## as.factor(VID)45983       0.437725   0.646546   0.677 0.498392    
## as.factor(VID)45984       0.687538   0.710773   0.967 0.333389    
## as.factor(VID)45985      -0.093996   0.576223  -0.163 0.870420    
## as.factor(VID)45986       1.169299   0.871620   1.342 0.179750    
## as.factor(VID)45988       0.909428   0.766977   1.186 0.235729    
## as.factor(VID)45992       1.994055   1.122148   1.777 0.075569 .  
## as.factor(VID)45994       0.700945   0.880713   0.796 0.426100    
## as.factor(VID)45997       0.243534   0.678435   0.359 0.719621    
## as.factor(VID)45999      15.794149 570.653349   0.028 0.977920    
## as.factor(VID)46000      15.754288 502.006275   0.031 0.974964    
## as.factor(VID)46003       0.318026   0.627836   0.507 0.612476    
## as.factor(VID)46004      -0.306279   0.608981  -0.503 0.615009    
## as.factor(VID)46006       0.043849   0.602844   0.073 0.942015    
## as.factor(VID)46007       0.865371   0.708115   1.222 0.221679    
## as.factor(VID)46008       0.007788   0.582599   0.013 0.989334    
## as.factor(VID)46010       2.367311   1.118584   2.116 0.034315 *  
## as.factor(VID)46076       0.756328   0.770048   0.982 0.326010    
## as.factor(VID)46079       0.342862   0.676593   0.507 0.612332    
## as.factor(VID)46084       0.272343   0.648838   0.420 0.674675    
## as.factor(VID)46091       0.173525   0.650558   0.267 0.789675    
## as.factor(VID)46092       0.081687   0.601186   0.136 0.891920    
## as.factor(VID)46095       1.658518   1.125874   1.473 0.140726    
## as.factor(VID)46100       0.838725   0.708661   1.184 0.236597    
## as.factor(VID)46109      -0.641287   0.572317  -1.121 0.262497    
## as.factor(VID)46110       1.213238   1.134402   1.069 0.284846    
## as.factor(VID)46113      -0.391143   0.641310  -0.610 0.541919    
## as.factor(VID)46116       1.435064   0.868034   1.653 0.098283 .  
## as.factor(VID)46117      -0.984451   0.570102  -1.727 0.084204 .  
## as.factor(VID)46120      -0.388318   0.589148  -0.659 0.509820    
## as.factor(VID)46127       0.065426   0.653773   0.100 0.920286    
## as.factor(VID)46128      -0.113219   0.634281  -0.178 0.858330    
## as.factor(VID)46135      -0.077296   0.616940  -0.125 0.900294    
## as.factor(VID)46138      -0.044558   0.655845  -0.068 0.945833    
## as.factor(VID)46139       1.369524   0.868723   1.576 0.114915    
## as.factor(VID)46145       0.331413   0.648231   0.511 0.609171    
## as.factor(VID)46146      -1.410913   0.540334  -2.611 0.009023 ** 
## as.factor(VID)46150       0.097166   0.631659   0.154 0.877747    
## as.factor(VID)46156      -0.076949   0.656446  -0.117 0.906685    
## as.factor(VID)46160       1.390067   0.868245   1.601 0.109375    
## as.factor(VID)46162       1.519163   1.127379   1.348 0.177814    
## as.factor(VID)46165       1.953333   1.121753   1.741 0.081627 .  
## as.factor(VID)46174      -0.019338   0.655750  -0.029 0.976474    
## as.factor(VID)46175      -0.092581   0.685307  -0.135 0.892538    
## as.factor(VID)46180       0.470813   0.713517   0.660 0.509351    
## as.factor(VID)46185       0.829648   0.769204   1.079 0.280775    
## as.factor(VID)46188       0.822161   0.768931   1.069 0.284968    
## as.factor(VID)46191      -0.964354   0.539430  -1.788 0.073820 .  
## as.factor(VID)46192      -0.840086   0.560803  -1.498 0.134132    
## as.factor(VID)46196      -0.901363   0.584039  -1.543 0.122751    
## as.factor(VID)46198      -0.583693   0.565771  -1.032 0.302224    
## as.factor(VID)46199      -0.246945   0.565983  -0.436 0.662611    
## as.factor(VID)46200      -0.869681   0.538006  -1.616 0.105989    
## as.factor(VID)46201       0.636847   0.671897   0.948 0.343214    
## as.factor(VID)46202      -0.701338   0.554632  -1.265 0.206047    
## as.factor(VID)46203       1.385665   0.762190   1.818 0.069064 .  
## as.factor(VID)46204       0.713032   0.670571   1.063 0.287637    
## as.factor(VID)46207      -1.191439   0.560819  -2.124 0.033631 *  
## as.factor(VID)46212       0.951454   0.706312   1.347 0.177957    
## as.factor(VID)46215      -0.853778   0.565914  -1.509 0.131383    
## as.factor(VID)46218      -0.007874   0.602813  -0.013 0.989579    
## as.factor(VID)46219       0.097076   0.590566   0.164 0.869434    
## as.factor(VID)46220      -0.380857   0.562571  -0.677 0.498410    
## as.factor(VID)46223       0.052711   0.602038   0.088 0.930230    
## as.factor(VID)46225       1.080990   0.705569   1.532 0.125502    
## as.factor(VID)46227      -0.649335   0.571987  -1.135 0.256280    
## as.factor(VID)46230      -0.470592   0.569454  -0.826 0.408582    
## as.factor(VID)46233      -0.084430   0.576263  -0.147 0.883516    
## as.factor(VID)46234       0.341177   0.626964   0.544 0.586323    
## as.factor(VID)46237      -0.940274   0.559284  -1.681 0.092722 .  
## as.factor(VID)46238       0.925564   0.668453   1.385 0.166164    
## as.factor(VID)46239      -1.275422   0.548323  -2.326 0.020016 *  
## as.factor(VID)46241      -0.308419   0.566757  -0.544 0.586316    
## as.factor(VID)46243       0.277143   0.628535   0.441 0.659260    
## as.factor(VID)46247       0.315299   0.627870   0.502 0.615546    
## as.factor(VID)46249       0.360611   0.676115   0.533 0.593786    
## as.factor(VID)46253      -0.008401   0.603372  -0.014 0.988891    
## as.factor(VID)46255       0.800840   0.669983   1.195 0.231964    
## as.factor(VID)46256       0.919532   0.707354   1.300 0.193615    
## as.factor(VID)46260       2.403916   1.117647   2.151 0.031486 *  
## as.factor(VID)46264      -0.463797   0.554987  -0.836 0.403330    
## as.factor(VID)46267       0.280572   0.627601   0.447 0.654835    
## as.factor(VID)46269      -0.450105   0.574969  -0.783 0.433725    
## as.factor(VID)46271      -0.523003   0.560011  -0.934 0.350347    
## as.factor(VID)46274      -0.005349   0.583025  -0.009 0.992680    
## as.factor(VID)46277      -0.810591   0.570047  -1.422 0.155035    
## as.factor(VID)46279       0.238918   0.612360   0.390 0.696418    
## as.factor(VID)46280       0.271716   0.628575   0.432 0.665543    
## as.factor(VID)46285       0.854080   0.642182   1.330 0.183530    
## as.factor(VID)46288       0.062731   0.582036   0.108 0.914172    
## as.factor(VID)46289       2.405344   1.117401   2.153 0.031348 *  
## as.factor(VID)46291       1.598472   0.866193   1.845 0.064979 .  
## as.factor(VID)46292      -0.298168   0.566528  -0.526 0.598675    
## as.factor(VID)46293      -0.413226   0.558307  -0.740 0.459214    
## as.factor(VID)46294       0.588207   0.672480   0.875 0.381747    
## as.factor(VID)46296      -0.300069   0.566833  -0.529 0.596543    
## as.factor(VID)46297       2.543307   1.116572   2.278 0.022740 *  
## as.factor(VID)46300       1.403154   0.762243   1.841 0.065647 .  
## as.factor(VID)46301      -0.130856   0.570239  -0.229 0.818499    
## as.factor(VID)46303       0.561799   0.644937   0.871 0.383704    
## as.factor(VID)46307      -1.047177   0.561393  -1.865 0.062137 .  
## as.factor(VID)46308      -0.670488   0.547399  -1.225 0.220627    
## as.factor(VID)46311       0.288017   0.628309   0.458 0.646664    
## as.factor(VID)46313      -0.255399   0.606883  -0.421 0.673874    
## as.factor(VID)46314       1.090252   0.764960   1.425 0.154088    
## as.factor(VID)46315       1.425087   0.867864   1.642 0.100577    
## as.factor(VID)46318       0.763148   0.708856   1.077 0.281663    
## as.factor(VID)46323      -0.474163   0.582336  -0.814 0.415505    
## as.factor(VID)46324      -0.484636   0.575678  -0.842 0.399870    
## as.factor(VID)46326       0.591139   0.644977   0.917 0.359391    
## as.factor(VID)46327       1.033604   0.705869   1.464 0.143112    
## as.factor(VID)46329       0.673635   0.671028   1.004 0.315434    
## as.factor(VID)46331      -0.249091   0.578909  -0.430 0.666994    
## as.factor(VID)46333       0.097083   0.631284   0.154 0.877778    
## as.factor(VID)46335       1.616169   0.865961   1.866 0.061995 .  
## as.factor(VID)46336       1.392701   0.762029   1.828 0.067606 .  
## as.factor(VID)46337       0.276409   0.612115   0.452 0.651584    
## as.factor(VID)46339      -0.769783   0.563699  -1.366 0.172068    
## as.factor(VID)46340      -0.334879   0.573030  -0.584 0.558951    
## as.factor(VID)46341       0.303506   0.648887   0.468 0.639975    
## as.factor(VID)46345      -0.521454   0.564290  -0.924 0.355440    
## as.factor(VID)46346      -0.022326   0.633389  -0.035 0.971881    
## as.factor(VID)46347      -0.962782   0.564268  -1.706 0.087961 .  
## as.factor(VID)46349      -0.333843   0.572639  -0.583 0.559900    
## as.factor(VID)46351       0.138599   0.630599   0.220 0.826035    
## as.factor(VID)46352      -0.474376   0.590468  -0.803 0.421749    
## as.factor(VID)46359       0.800085   0.669724   1.195 0.232224    
## as.factor(VID)46363       0.166418   0.600371   0.277 0.781633    
## as.factor(VID)46364       0.896745   0.767713   1.168 0.242777    
## as.factor(VID)46367      -0.274640   0.595970  -0.461 0.644922    
## as.factor(VID)46371       0.199001   0.612696   0.325 0.745336    
## as.factor(VID)46375       0.885464   0.707895   1.251 0.210993    
## as.factor(VID)46377       0.322650   0.598707   0.539 0.589948    
## as.factor(VID)46378       1.633888   0.866294   1.886 0.059286 .  
## as.factor(VID)46381       0.477345   0.597337   0.799 0.424221    
## as.factor(VID)46382      -0.758383   0.552227  -1.373 0.169653    
## as.factor(VID)46383      15.753441 506.080935   0.031 0.975167    
## as.factor(VID)46390       2.357532   1.118022   2.109 0.034974 *  
## as.factor(VID)46392      -0.336865   0.609882  -0.552 0.580713    
## as.factor(VID)46394       1.587244   0.866414   1.832 0.066956 .  
## as.factor(VID)46399       1.568938   0.866335   1.811 0.070140 .  
## as.factor(VID)46416      -0.631737   0.578057  -1.093 0.274454    
## as.factor(VID)46419      -0.606112   0.565855  -1.071 0.284105    
## as.factor(VID)46420       0.685458   0.710447   0.965 0.334631    
## as.factor(VID)46423      -0.618386   0.566242  -1.092 0.274795    
## as.factor(VID)46424       0.719185   0.670510   1.073 0.283453    
## as.factor(VID)46426      -0.649875   0.557848  -1.165 0.244032    
## as.factor(VID)46428       1.635455   0.865935   1.889 0.058938 .  
## as.factor(VID)46430      -0.348572   0.572746  -0.609 0.542791    
## as.factor(VID)46431       0.096297   0.651960   0.148 0.882577    
## as.factor(VID)46436      -0.406665   0.581403  -0.699 0.484268    
## as.factor(VID)46437      -0.336485   0.567036  -0.593 0.552906    
## as.factor(VID)46440      15.727451 548.287558   0.029 0.977116    
## as.factor(VID)46442       0.543956   0.645151   0.843 0.399147    
## as.factor(VID)46464       0.858625   0.770823   1.114 0.265319    
## as.factor(VID)46465       0.019320   0.615335   0.031 0.974952    
## as.factor(VID)46483      -0.550780   0.565595  -0.974 0.330152    
## as.factor(VID)46485      -0.106627   0.584074  -0.183 0.855146    
## as.factor(VID)46492      -0.275214   0.596396  -0.461 0.644467    
## as.factor(VID)46514       0.256792   0.679878   0.378 0.705651    
## as.factor(VID)46545       0.602108   0.672373   0.895 0.370522    
## as.factor(VID)46569       0.165345   0.680159   0.243 0.807930    
## as.factor(VID)46583       0.280108   0.649502   0.431 0.666274    
## as.factor(VID)46592       1.110269   0.765124   1.451 0.146753    
## as.factor(VID)46598       0.168815   0.651541   0.259 0.795557    
## as.factor(VID)46600       0.506586   0.712593   0.711 0.477143    
## as.factor(VID)46607      15.789187 532.961827   0.030 0.976366    
## as.factor(VID)46612      -0.707024   0.554269  -1.276 0.202098    
## as.factor(VID)46613      -0.582598   0.593334  -0.982 0.326146    
## as.factor(VID)46615      15.856998 668.000943   0.024 0.981062    
## as.factor(VID)46625       1.373681   0.868538   1.582 0.113741    
## as.factor(VID)46628      -0.654854   0.594678  -1.101 0.270814    
## as.factor(VID)46643      -0.022523   0.654512  -0.034 0.972548    
## as.factor(VID)46670       0.352218   0.716407   0.492 0.622970    
## as.factor(VID)46779       1.437921   0.761894   1.887 0.059120 .  
## as.factor(VID)46810       0.364714   0.647937   0.563 0.573513    
## as.factor(VID)46818       0.145724   0.630605   0.231 0.817248    
## as.factor(VID)46823      -0.070053   0.583372  -0.120 0.904417    
## as.factor(VID)46852       0.208313   0.614285   0.339 0.734523    
## as.factor(VID)46859       0.464362   0.626758   0.741 0.458757    
## as.factor(VID)46866       2.152194   1.119490   1.922 0.054546 .  
## as.factor(VID)46875       0.147566   0.590636   0.250 0.802710    
## as.factor(VID)46904       0.077549   0.654037   0.119 0.905616    
## as.factor(VID)46908      -0.353299   0.580058  -0.609 0.542474    
## as.factor(VID)46912      -1.137630   0.578899  -1.965 0.049396 *  
## as.factor(VID)46917       0.860185   0.669340   1.285 0.198749    
## as.factor(VID)46919       0.632666   0.711825   0.889 0.374113    
## as.factor(VID)46924      -0.308888   0.598228  -0.516 0.605618    
## as.factor(VID)46925       1.143535   0.705239   1.621 0.104914    
## as.factor(VID)46928       0.101189   0.613968   0.165 0.869092    
## as.factor(VID)46931      -0.804163   0.560953  -1.434 0.151696    
## as.factor(VID)46932      -0.310913   0.579439  -0.537 0.591561    
## as.factor(VID)46933       2.085552   1.120168   1.862 0.062628 .  
## as.factor(VID)46941      -0.908366   0.558083  -1.628 0.103598    
## as.factor(VID)46942      -0.808737   0.599481  -1.349 0.177317    
## as.factor(VID)46943       0.110100   0.680759   0.162 0.871517    
## as.factor(VID)46944       1.298143   0.762875   1.702 0.088822 .  
## as.factor(VID)46946       0.572971   0.625077   0.917 0.359331    
## as.factor(VID)46956      15.780458 558.986502   0.028 0.977478    
## as.factor(VID)46958       0.633055   0.711878   0.889 0.373856    
## as.factor(VID)46959       0.209987   0.600073   0.350 0.726386    
## as.factor(VID)46961      -1.132175   0.554174  -2.043 0.041053 *  
## as.factor(VID)46962       1.060847   0.765490   1.386 0.165796    
## as.factor(VID)46966       0.532160   0.625005   0.851 0.394520    
## as.factor(VID)46968       1.518539   0.867394   1.751 0.079999 .  
## as.factor(VID)46970       1.595960   0.866425   1.842 0.065474 .  
## as.factor(VID)46973       0.418372   0.626510   0.668 0.504273    
## as.factor(VID)46975       1.605362   0.866887   1.852 0.064045 .  
## as.factor(VID)46976       0.121130   0.590994   0.205 0.837603    
## as.factor(VID)46981       0.628149   0.644642   0.974 0.329851    
## as.factor(VID)46984       0.467598   0.626099   0.747 0.455158    
## as.factor(VID)46986      -0.195204   0.585239  -0.334 0.738722    
## as.factor(VID)46988       0.240343   0.612833   0.392 0.694923    
## as.factor(VID)46991       0.202065   0.589864   0.343 0.731928    
## as.factor(VID)47025       2.232469   1.118523   1.996 0.045944 *  
## as.factor(VID)47041      -0.704012   0.548189  -1.284 0.199054    
## as.factor(VID)47043      -0.413714   0.626311  -0.661 0.508897    
## as.factor(VID)47108       0.276550   0.628373   0.440 0.659861    
## as.factor(VID)47120       1.990463   1.121564   1.775 0.075944 .  
## as.factor(VID)47145      -0.244705   0.561061  -0.436 0.662731    
## as.factor(VID)47146      -0.456663   0.551549  -0.828 0.407691    
## as.factor(VID)47150      -0.598387   0.552671  -1.083 0.278934    
## as.factor(VID)47153      -0.390709   0.558023  -0.700 0.483824    
## as.factor(VID)47154       0.186194   0.600152   0.310 0.756375    
## as.factor(VID)47156       0.602640   0.644610   0.935 0.349844    
## as.factor(VID)47157      -0.370700   0.557864  -0.664 0.506370    
## as.factor(VID)47158       0.099284   0.582261   0.171 0.864606    
## as.factor(VID)47160       0.192511   0.599937   0.321 0.748297    
## as.factor(VID)47291      -0.153570   0.593777  -0.259 0.795918    
## as.factor(VID)47295      -1.102737   0.567712  -1.942 0.052086 .  
## as.factor(VID)47296       2.487725   1.116943   2.227 0.025930 *  
## as.factor(VID)47297      -0.583901   0.577203  -1.012 0.311728    
## as.factor(VID)47298       0.386300   0.676141   0.571 0.567776    
## as.factor(VID)47299      -0.085749   0.576048  -0.149 0.881666    
## as.factor(VID)47300      -0.964665   0.587783  -1.641 0.100757    
## as.factor(VID)47301       0.088483   0.601117   0.147 0.882975    
## as.factor(VID)47302       2.257053   1.118955   2.017 0.043684 *  
## as.factor(VID)47303      -1.764014   0.523868  -3.367 0.000759 ***
## as.factor(VID)47304      15.789555 518.895611   0.030 0.975725    
## as.factor(VID)47305       0.283663   0.648692   0.437 0.661905    
## as.factor(VID)47306      -0.163797   0.595391  -0.275 0.783233    
## as.factor(VID)47308       0.778488   0.669757   1.162 0.245096    
## as.factor(VID)47312      -0.805214   0.560241  -1.437 0.150643    
## as.factor(VID)47313       0.955393   0.706826   1.352 0.176482    
## as.factor(VID)47318      -0.345227   0.562451  -0.614 0.539354    
## as.factor(VID)47327       2.562371   1.116831   2.294 0.021772 *  
## as.factor(VID)47332       0.392205   0.611203   0.642 0.521073    
## as.factor(VID)47338       0.070826   0.601863   0.118 0.906323    
## as.factor(VID)47339      -1.521528   0.538599  -2.825 0.004728 ** 
## as.factor(VID)47345      15.805475 519.074198   0.030 0.975709    
## as.factor(VID)47352       0.149910   0.581358   0.258 0.796514    
## as.factor(VID)47355      15.722019 486.682243   0.032 0.974229    
## as.factor(VID)47356       0.492063   0.646052   0.762 0.446271    
## as.factor(VID)47409      -0.531779   0.570153  -0.933 0.350977    
## as.factor(VID)47582      -0.658947   0.551248  -1.195 0.231941    
## as.factor(VID)57619      -0.072776   0.604272  -0.120 0.904138    
## as.factor(VID)57636      15.796454 523.532430   0.030 0.975929    
## as.factor(VID)57641       2.176013   1.119156   1.944 0.051855 .  
## as.factor(VID)57643       2.295024   1.118420   2.052 0.040167 *  
## as.factor(VID)57644       0.568646   0.672774   0.845 0.397985    
## as.factor(VID)57646      -0.295638   0.597261  -0.495 0.620607    
## as.factor(VID)57648       0.837509   0.669554   1.251 0.210991    
## as.factor(VID)67891       0.031481   0.632009   0.050 0.960273    
## as.factor(VID)67936      15.775771 641.018933   0.025 0.980366    
## as.factor(VID)67944      -0.125314   0.618708  -0.203 0.839494    
## as.factor(VID)67952       1.560857   1.126862   1.385 0.166011    
## as.factor(VID)67958      15.781843 582.770380   0.027 0.978395    
## as.factor(VID)67961      15.821342 678.043819   0.023 0.981384    
## as.factor(VID)67963      -0.318569   0.597643  -0.533 0.594004    
## as.factor(VID)67965       0.977977   0.766847   1.275 0.202195    
## as.factor(VID)68080       1.127498   0.764962   1.474 0.140501    
## as.factor(VID)68083       1.081570   0.765144   1.414 0.157493    
## as.factor(VID)78079       1.374272   0.762279   1.803 0.071412 .  
## as.factor(VID)78205      15.784982 589.157368   0.027 0.978625    
## as.factor(VID)78214       0.026923   0.632497   0.043 0.966048    
## as.factor(VID)78215       1.550169   0.866733   1.789 0.073692 .  
## as.factor(VID)78231       2.288708   1.118076   2.047 0.040657 *  
## as.factor(VID)78232      -0.372018   0.573657  -0.649 0.516659    
## as.factor(VID)78233      15.739795 498.063149   0.032 0.974789    
## as.factor(VID)78234       0.301659   0.629624   0.479 0.631861    
## as.factor(VID)78235      15.744419 506.152991   0.031 0.975185    
## as.factor(VID)78236      -0.709223   0.563773  -1.258 0.208394    
## as.factor(VID)78237       0.004801   0.591431   0.008 0.993523    
## as.factor(VID)78238      15.752656 514.709173   0.031 0.975585    
## as.factor(VID)78239      -1.988252   0.550100  -3.614 0.000301 ***
## as.factor(VID)78240      -1.310650   0.556127  -2.357 0.018436 *  
## as.factor(VID)78241      15.743716 510.359310   0.031 0.975391    
## as.factor(VID)78242       0.872311   0.709179   1.230 0.218686    
## as.factor(VID)78243      -1.283590   0.530650  -2.419 0.015567 *  
## as.factor(VID)78244      -1.246037   0.536366  -2.323 0.020173 *  
## as.factor(VID)78245      15.743881 523.589969   0.030 0.976012    
## as.factor(VID)78246      -1.106187   0.549741  -2.012 0.044199 *  
## as.factor(VID)78247      15.708400 416.447914   0.038 0.969911    
## as.factor(VID)78248      -0.033518   0.602772  -0.056 0.955655    
## as.factor(VID)78249      -0.995070   0.547815  -1.816 0.069304 .  
## as.factor(VID)78250      15.785718 486.524569   0.032 0.974116    
## as.factor(VID)78251      -1.576847   0.561215  -2.810 0.004959 ** 
## as.factor(VID)78252      -0.695917   0.581068  -1.198 0.231052    
## as.factor(VID)78253       2.262027   1.118899   2.022 0.043212 *  
## as.factor(VID)78254       1.751807   0.865985   2.023 0.043083 *  
## as.factor(VID)78255      15.830001 589.309921   0.027 0.978570    
## as.factor(VID)78256      15.863330 564.771063   0.028 0.977592    
## as.factor(VID)78257      15.791484 459.523260   0.034 0.972586    
## as.factor(VID)78258      -2.267406   0.547946  -4.138 3.50e-05 ***
## as.factor(VID)78259      -1.129091   0.563852  -2.002 0.045235 *  
## as.factor(VID)78260      -0.243663   0.596842  -0.408 0.683088    
## as.factor(VID)78261      -0.293771   0.623220  -0.471 0.637373    
## as.factor(VID)78262      15.818045 688.251716   0.023 0.981664    
## as.factor(VID)78263       2.338700   1.117720   2.092 0.036404 *  
## as.factor(VID)78266       0.087066   0.602409   0.145 0.885083    
## as.factor(VID)78267       0.842959   0.769826   1.095 0.273517    
## as.factor(VID)78268      -1.412754   0.545757  -2.589 0.009636 ** 
## as.factor(VID)78269      -1.541264   0.542124  -2.843 0.004469 ** 
## as.factor(VID)78270       2.481996   1.117535   2.221 0.026354 *  
## as.factor(VID)78271      15.821523 486.556290   0.033 0.974059    
## as.factor(VID)78272       2.299845   1.119029   2.055 0.039858 *  
## as.factor(VID)78273       1.385150   0.762101   1.818 0.069134 .  
## as.factor(VID)78276       1.780017   0.865352   2.057 0.039688 *  
## as.factor(VID)78277      -0.338993   0.598968  -0.566 0.571420    
## as.factor(VID)78279       0.273566   0.599968   0.456 0.648413    
## as.factor(VID)78283       0.214436   0.589283   0.364 0.715938    
## as.factor(VID)78285      15.803299 490.303674   0.032 0.974287    
## as.factor(VID)78286      15.791458 447.561366   0.035 0.971854    
## as.factor(VID)78288       2.346093   1.118685   2.097 0.035977 *  
## as.factor(VID)78291       0.140548   0.601302   0.234 0.815188    
## as.factor(VID)78292      15.795047 625.050119   0.025 0.979840    
## as.factor(VID)78295       2.389088   1.118482   2.136 0.032679 *  
## as.factor(VID)78298      15.761151 510.187737   0.031 0.975355    
## as.factor(VID)78300       1.620474   0.866619   1.870 0.061500 .  
## as.factor(VID)78301       1.554667   0.866750   1.794 0.072865 .  
## as.factor(VID)78302      -0.927305   0.558781  -1.660 0.097012 .  
## as.factor(VID)78303       2.248605   1.119261   2.009 0.044536 *  
## as.factor(VID)78309       0.186767   0.722130   0.259 0.795918    
## as.factor(VID)78312       0.642393   0.671884   0.956 0.339018    
## as.factor(VID)78313      15.753804 514.346726   0.031 0.975566    
## as.factor(VID)78542      -0.711497   0.575134  -1.237 0.216051    
## as.factor(VID)78592       0.892926   0.707791   1.262 0.207105    
## as.factor(VID)78593       0.258996   0.599295   0.432 0.665619    
## as.factor(VID)78594       1.655915   0.866624   1.911 0.056035 .  
## as.factor(VID)78595       0.600128   0.645188   0.930 0.352289    
## as.factor(VID)78596      -0.057576   0.584194  -0.099 0.921491    
## as.factor(VID)78597       0.272549   0.628454   0.434 0.664520    
## as.factor(VID)78598       1.237141   0.764049   1.619 0.105406    
## as.factor(VID)78599       0.992409   0.707280   1.403 0.160576    
## as.factor(VID)78600      -0.119587   0.659266  -0.181 0.856058    
## as.factor(VID)78601       0.618837   0.624635   0.991 0.321823    
## as.factor(VID)78602      -0.117786   0.619920  -0.190 0.849307    
## as.factor(VID)78604      -0.026189   0.633961  -0.041 0.967049    
## as.factor(VID)78605      -1.527879   0.537106  -2.845 0.004446 ** 
## as.factor(VID)78606       0.260018   0.599983   0.433 0.664743    
## as.factor(VID)78607       1.617264   0.866599   1.866 0.062011 .  
## as.factor(VID)78608      -1.831664   0.578401  -3.167 0.001541 ** 
## as.factor(VID)78609       0.037324   0.582709   0.064 0.948929    
## as.factor(VID)78611      -0.132736   0.570417  -0.233 0.815994    
## as.factor(VID)78613       0.032167   0.603057   0.053 0.957461    
## as.factor(VID)78614       1.646500   1.126585   1.461 0.143879    
## as.factor(VID)78616      -1.151708   0.588214  -1.958 0.050233 .  
## as.factor(VID)78617       0.808014   0.710055   1.138 0.255137    
## as.factor(VID)78619       0.939569   0.767507   1.224 0.220883    
## as.factor(VID)78638       0.743653   0.710549   1.047 0.295289    
## as.factor(VID)78643       0.084692   0.581771   0.146 0.884256    
## as.factor(VID)78649      -1.113646   0.546270  -2.039 0.041486 *  
## as.factor(VID)78650       1.422021   0.868324   1.638 0.101492    
## as.factor(VID)78651       1.302766   0.703869   1.851 0.064189 .  
## as.factor(VID)78652      -0.598458   0.566787  -1.056 0.291024    
## as.factor(VID)78654      -0.366908   0.609070  -0.602 0.546904    
## as.factor(VID)78683       2.355678   1.117978   2.107 0.035110 *  
## as.factor(VID)78687       0.833156   0.708954   1.175 0.239919    
## as.factor(VID)78689       0.168448   0.613980   0.274 0.783812    
## as.factor(VID)78691       0.708219   0.671415   1.055 0.291509    
## as.factor(VID)78693       0.580813   0.646215   0.899 0.368764    
## as.factor(VID)78694      -0.107306   0.592706  -0.181 0.856332    
## as.factor(VID)78695       0.085388   0.614725   0.139 0.889526    
## as.factor(VID)78696       0.225318   0.629479   0.358 0.720385    
## as.factor(VID)78698       1.345215   0.762844   1.763 0.077829 .  
## as.factor(VID)78699       0.021065   0.615911   0.034 0.972717    
## as.factor(VID)78700       0.382819   0.676908   0.566 0.571706    
## as.factor(VID)78701      -0.482984   0.583569  -0.828 0.407875    
## as.factor(VID)78703       0.555820   0.645901   0.861 0.389494    
## as.factor(VID)78704       1.114903   0.705257   1.581 0.113913    
## as.factor(VID)78705       0.442459   0.625850   0.707 0.479584    
## as.factor(VID)78706       0.363225   0.648378   0.560 0.575339    
## as.factor(VID)78708       0.059886   0.632239   0.095 0.924536    
## as.factor(VID)78709       0.166935   0.651516   0.256 0.797777    
## as.factor(VID)78711       0.703494   0.670946   1.049 0.294403    
## as.factor(VID)78712       0.108207   0.630624   0.172 0.863762    
## as.factor(VID)78713       1.321680   0.869359   1.520 0.128437    
## as.factor(VID)78717       1.671492   0.865589   1.931 0.053477 .  
## as.factor(VID)78718       0.189608   0.629416   0.301 0.763227    
## as.factor(VID)78739       0.275108   0.628039   0.438 0.661355    
## as.factor(VID)78741       0.452720   0.647100   0.700 0.484169    
## as.factor(VID)78742      -0.344545   0.588934  -0.585 0.558526    
## as.factor(VID)78744      -0.066329   0.604724  -0.110 0.912660    
## as.factor(VID)78745       0.187966   0.614178   0.306 0.759570    
## as.factor(VID)78746       0.215480   0.629842   0.342 0.732263    
## as.factor(VID)78748       0.116072   0.614800   0.189 0.850253    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 37106  on 47598  degrees of freedom
## Residual deviance: 32845  on 46789  degrees of freedom
## AIC: 34465
## 
## Number of Fisher Scoring iterations: 16
## Warning in ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary =
## TRUE, : maximum three independent variables are allowed
## NULL
Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 6.8991, df = 8, p-value = 0.5476
148
glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 148))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01, family = "binomial", data = child_ica_dummy %>% 
##     filter(DID == 148))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.1602   0.4780   0.5202   0.5651   0.6721  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              2.08964    0.22845   9.147   <2e-16 ***
## n_children_in_household -0.05985    0.04872  -1.228    0.219    
## PR004_PR009_01           0.20142    0.15237   1.322    0.186    
## C002_01                 -0.29794    0.13730  -2.170    0.030 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1462.2  on 1851  degrees of freedom
## Residual deviance: 1454.1  on 1848  degrees of freedom
## AIC: 1462.1
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 0.70355, df = 8, p-value = 0.9995
151
glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 151))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01, family = "binomial", data = child_ica_dummy %>% 
##     filter(DID == 151))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.0776   0.5003   0.5051   0.5338   0.7164  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              1.42858    0.25219   5.665 1.47e-08 ***
## n_children_in_household -0.01015    0.05416  -0.187 0.851379    
## PR004_PR009_01           0.61702    0.16072   3.839 0.000123 ***
## C002_01                 -0.12832    0.13800  -0.930 0.352425    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1445.1  on 1793  degrees of freedom
## Residual deviance: 1430.1  on 1790  degrees of freedom
## AIC: 1438.1
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 26.905, df = 8, p-value = 0.0007342
260
glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01, family = "binomial", data = child_ica_dummy %>% filter(DID == 260))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01, family = "binomial", data = child_ica_dummy %>% 
##     filter(DID == 260))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.1956   0.4753   0.5348   0.5847   0.7639  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              2.48695    0.21355  11.646   <2e-16 ***
## n_children_in_household -0.08535    0.03699  -2.308   0.0210 *  
## PR004_PR009_01          -0.27317    0.14617  -1.869   0.0616 .  
## C002_01                 -0.27779    0.12815  -2.168   0.0302 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1650.2  on 2004  degrees of freedom
## Residual deviance: 1636.2  on 2001  degrees of freedom
## AIC: 1644.2
## 
## Number of Fisher Scoring iterations: 4
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

## Warning in eval(family$initialize): non-integer #successes in a binomial glm!

Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 15.649, df = 8, p-value = 0.0477
260, VID
glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + factor(VID), family = "binomial", data = child_ica_dummy %>% filter(DID == 260))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01 + factor(VID), family = "binomial", data = child_ica_dummy %>% 
##     filter(DID == 260))
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.98580   0.00014   0.42879   0.57612   1.18235  
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               1.67422    0.48297   3.466 0.000527 ***
## n_children_in_household   0.03091    0.04526   0.683 0.494630    
## PR004_PR009_01            0.06956    0.18637   0.373 0.708963    
## C002_01                  -0.12462    0.14000  -0.890 0.373392    
## factor(VID)47295         -1.07596    0.53590  -2.008 0.044668 *  
## factor(VID)47296          2.51948    1.10231   2.286 0.022276 *  
## factor(VID)47297         -0.54884    0.54600  -1.005 0.314801    
## factor(VID)47298          0.40838    0.64832   0.630 0.528759    
## factor(VID)47299         -0.03990    0.55358  -0.072 0.942536    
## factor(VID)47300         -0.94556    0.55578  -1.701 0.088881 .  
## factor(VID)47301          0.12755    0.57501   0.222 0.824453    
## factor(VID)47302          2.28506    1.10787   2.063 0.039154 *  
## factor(VID)47303         -1.72351    0.50253  -3.430 0.000604 ***
## factor(VID)47304         16.81541  855.93786   0.020 0.984326    
## factor(VID)47305          0.31982    0.62245   0.514 0.607388    
## factor(VID)47306         -0.13887    0.57337  -0.242 0.808632    
## factor(VID)47308          0.81913    0.64618   1.268 0.204919    
## factor(VID)47312         -0.76830    0.53045  -1.448 0.147507    
## factor(VID)47313          0.98803    0.68135   1.450 0.147027    
## factor(VID)47318         -0.31342    0.53646  -0.584 0.559061    
## factor(VID)47327          2.58615    1.10413   2.342 0.019168 *  
## factor(VID)47332          0.42155    0.59249   0.711 0.476779    
## factor(VID)47338          0.09982    0.57314   0.174 0.861735    
## factor(VID)47339         -1.50129    0.51737  -2.902 0.003711 ** 
## factor(VID)47345         16.82728  856.13930   0.020 0.984319    
## factor(VID)47352          0.19564    0.56564   0.346 0.729437    
## factor(VID)47355         16.76490  802.63292   0.021 0.983335    
## factor(VID)47356          0.52232    0.62125   0.841 0.400484    
## factor(VID)47409         -0.49782    0.53881  -0.924 0.355524    
## factor(VID)47582         -0.63324    0.52455  -1.207 0.227346    
## factor(VID)57636         16.81954  863.54271   0.019 0.984460    
## factor(VID)57643          2.31954    1.10185   2.105 0.035279 *  
## factor(VID)57644          0.59762    0.64630   0.925 0.355130    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1650.2  on 2004  degrees of freedom
## Residual deviance: 1365.7  on 1972  degrees of freedom
## AIC: 1431.7
## 
## Number of Fisher Scoring iterations: 17
## Warning in ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary =
## TRUE, : maximum three independent variables are allowed
## NULL
Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 2.569, df = 8, p-value = 0.9584

Without dists_not_fit. GLM C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(DID)

glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01 + as.factor(DID), family = "binomial", data = child_ica_dummy %>% filter(!DID %in% dists_not_fit))

ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
## 
## Call:
## glm(formula = C003_01 ~ n_children_in_household + PR004_PR009_01 + 
##     C002_01 + as.factor(DID), family = "binomial", data = child_ica_dummy %>% 
##     filter(!DID %in% dists_not_fit))
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.5291  -1.1412   0.6592   0.8972   1.7474  
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)              2.20023    0.08886  24.761  < 2e-16 ***
## n_children_in_household -0.03221    0.00342  -9.417  < 2e-16 ***
## PR004_PR009_01           0.36964    0.01108  33.351  < 2e-16 ***
## C002_01                 -0.73250    0.01017 -72.036  < 2e-16 ***
## as.factor(DID)147       -0.61110    0.10692  -5.715 1.09e-08 ***
## as.factor(DID)149       -0.80774    0.10286  -7.853 4.07e-15 ***
## as.factor(DID)150       -0.78424    0.10201  -7.688 1.50e-14 ***
## as.factor(DID)153       -0.90444    0.10406  -8.692  < 2e-16 ***
## as.factor(DID)154       -0.68180    0.10568  -6.452 1.11e-10 ***
## as.factor(DID)155       -0.84481    0.10119  -8.349  < 2e-16 ***
## as.factor(DID)156       -1.33487    0.11582 -11.525  < 2e-16 ***
## as.factor(DID)157       -0.45598    0.11098  -4.108 3.98e-05 ***
## as.factor(DID)159       -1.32982    0.10655 -12.481  < 2e-16 ***
## as.factor(DID)160       -0.73750    0.10489  -7.031 2.05e-12 ***
## as.factor(DID)161       -1.27916    0.09762 -13.103  < 2e-16 ***
## as.factor(DID)164       -0.57774    0.11384  -5.075 3.87e-07 ***
## as.factor(DID)165       -0.98013    0.10429  -9.398  < 2e-16 ***
## as.factor(DID)166       -0.83499    0.10473  -7.973 1.55e-15 ***
## as.factor(DID)167        0.50880    0.13040   3.902 9.55e-05 ***
## as.factor(DID)169       -1.72476    0.09928 -17.373  < 2e-16 ***
## as.factor(DID)170       -0.47695    0.10872  -4.387 1.15e-05 ***
## as.factor(DID)171       -0.44109    0.11072  -3.984 6.78e-05 ***
## as.factor(DID)172       -0.63932    0.10654  -6.001 1.96e-09 ***
## as.factor(DID)174       -0.96366    0.10561  -9.125  < 2e-16 ***
## as.factor(DID)175       -1.33897    0.10151 -13.191  < 2e-16 ***
## as.factor(DID)176       -0.48230    0.11479  -4.201 2.65e-05 ***
## as.factor(DID)177       -0.85401    0.10314  -8.280  < 2e-16 ***
## as.factor(DID)179       -0.72202    0.10558  -6.839 7.98e-12 ***
## as.factor(DID)180       -0.68083    0.10615  -6.414 1.42e-10 ***
## as.factor(DID)182       -1.39805    0.10043 -13.920  < 2e-16 ***
## as.factor(DID)183       -0.65956    0.10778  -6.120 9.39e-10 ***
## as.factor(DID)184       -1.30777    0.10005 -13.071  < 2e-16 ***
## as.factor(DID)186       -1.69771    0.09953 -17.058  < 2e-16 ***
## as.factor(DID)187       -1.46210    0.10117 -14.452  < 2e-16 ***
## as.factor(DID)188       -1.30548    0.10037 -13.006  < 2e-16 ***
## as.factor(DID)189       -0.35890    0.10834  -3.313 0.000924 ***
## as.factor(DID)190       -0.99623    0.10865  -9.169  < 2e-16 ***
## as.factor(DID)191       -0.90900    0.10845  -8.382  < 2e-16 ***
## as.factor(DID)192       -1.06589    0.10269 -10.380  < 2e-16 ***
## as.factor(DID)194       -1.68630    0.09992 -16.876  < 2e-16 ***
## as.factor(DID)195       -1.88739    0.10107 -18.674  < 2e-16 ***
## as.factor(DID)196       -0.70672    0.10997  -6.426 1.31e-10 ***
## as.factor(DID)197       -1.43308    0.10084 -14.211  < 2e-16 ***
## as.factor(DID)198       -1.79500    0.09975 -17.996  < 2e-16 ***
## as.factor(DID)200       -1.46069    0.10294 -14.190  < 2e-16 ***
## as.factor(DID)202       -1.48817    0.10096 -14.740  < 2e-16 ***
## as.factor(DID)203       -1.44779    0.09921 -14.593  < 2e-16 ***
## as.factor(DID)204       -0.89171    0.10343  -8.621  < 2e-16 ***
## as.factor(DID)205       -1.59378    0.09841 -16.196  < 2e-16 ***
## as.factor(DID)206       -1.65528    0.09837 -16.827  < 2e-16 ***
## as.factor(DID)207       -2.02987    0.09995 -20.309  < 2e-16 ***
## as.factor(DID)208       -1.51565    0.09950 -15.233  < 2e-16 ***
## as.factor(DID)209       -1.03418    0.10817  -9.560  < 2e-16 ***
## as.factor(DID)210       -1.72023    0.09790 -17.572  < 2e-16 ***
## as.factor(DID)211       -1.03761    0.10795  -9.611  < 2e-16 ***
## as.factor(DID)212       -1.34488    0.10228 -13.149  < 2e-16 ***
## as.factor(DID)213       -1.58217    0.09868 -16.033  < 2e-16 ***
## as.factor(DID)214       -2.10360    0.09920 -21.205  < 2e-16 ***
## as.factor(DID)215       -1.35788    0.09697 -14.003  < 2e-16 ***
## as.factor(DID)216       -1.68497    0.10160 -16.585  < 2e-16 ***
## as.factor(DID)217       -1.71022    0.09997 -17.107  < 2e-16 ***
## as.factor(DID)218       -1.64404    0.09735 -16.889  < 2e-16 ***
## as.factor(DID)219       -1.54732    0.10069 -15.366  < 2e-16 ***
## as.factor(DID)220       -2.09984    0.09644 -21.774  < 2e-16 ***
## as.factor(DID)221       -2.45952    0.10142 -24.251  < 2e-16 ***
## as.factor(DID)222       -1.69695    0.09785 -17.341  < 2e-16 ***
## as.factor(DID)223       -1.61343    0.10554 -15.288  < 2e-16 ***
## as.factor(DID)224       -2.15683    0.09912 -21.759  < 2e-16 ***
## as.factor(DID)225       -1.86448    0.09974 -18.693  < 2e-16 ***
## as.factor(DID)226       -1.38124    0.10295 -13.417  < 2e-16 ***
## as.factor(DID)227       -1.67171    0.09806 -17.047  < 2e-16 ***
## as.factor(DID)228       -2.25942    0.10320 -21.893  < 2e-16 ***
## as.factor(DID)229       -1.96803    0.10194 -19.306  < 2e-16 ***
## as.factor(DID)230       -2.14701    0.09950 -21.579  < 2e-16 ***
## as.factor(DID)231       -1.40315    0.09889 -14.188  < 2e-16 ***
## as.factor(DID)233       -2.27841    0.10437 -21.830  < 2e-16 ***
## as.factor(DID)234       -1.58568    0.10526 -15.065  < 2e-16 ***
## as.factor(DID)235       -0.83483    0.10293  -8.111 5.02e-16 ***
## as.factor(DID)236       -1.10526    0.10314 -10.717  < 2e-16 ***
## as.factor(DID)237       -0.76287    0.10610  -7.190 6.46e-13 ***
## as.factor(DID)239       -0.89327    0.10257  -8.708  < 2e-16 ***
## as.factor(DID)240       -1.24717    0.10232 -12.189  < 2e-16 ***
## as.factor(DID)241       -2.37464    0.09837 -24.140  < 2e-16 ***
## as.factor(DID)242       -1.09923    0.10351 -10.619  < 2e-16 ***
## as.factor(DID)243       -1.07306    0.10870  -9.872  < 2e-16 ***
## as.factor(DID)244       -0.56597    0.11096  -5.101 3.38e-07 ***
## as.factor(DID)246       -1.71773    0.10021 -17.142  < 2e-16 ***
## as.factor(DID)247       -1.23790    0.10783 -11.480  < 2e-16 ***
## as.factor(DID)248       -0.94198    0.09993  -9.427  < 2e-16 ***
## as.factor(DID)250       -0.72232    0.10526  -6.862 6.77e-12 ***
## as.factor(DID)251       -0.68896    0.10687  -6.447 1.14e-10 ***
## as.factor(DID)252       -0.91843    0.10390  -8.840  < 2e-16 ***
## as.factor(DID)253       -0.51896    0.11209  -4.630 3.66e-06 ***
## as.factor(DID)254       -1.76683    0.11543 -15.306  < 2e-16 ***
## as.factor(DID)255       -0.86631    0.10906  -7.943 1.97e-15 ***
## as.factor(DID)257        0.49870    0.14055   3.548 0.000388 ***
## as.factor(DID)258       -1.57776    0.09915 -15.912  < 2e-16 ***
## as.factor(DID)259       -0.84587    0.10732  -7.882 3.23e-15 ***
## as.factor(DID)261       -2.30250    0.09833 -23.416  < 2e-16 ***
## as.factor(DID)262       -0.72491    0.10371  -6.990 2.75e-12 ***
## as.factor(DID)263       -0.41112    0.10755  -3.822 0.000132 ***
## as.factor(DID)264       -0.44120    0.10580  -4.170 3.04e-05 ***
## as.factor(DID)266        0.61879    0.13022   4.752 2.02e-06 ***
## as.factor(DID)268       -1.12253    0.10429 -10.763  < 2e-16 ***
## as.factor(DID)275       -1.11609    0.10208 -10.933  < 2e-16 ***
## as.factor(DID)278       -0.79219    0.09862  -8.033 9.50e-16 ***
## as.factor(DID)279       -1.83698    0.10016 -18.340  < 2e-16 ***
## as.factor(DID)280       -1.05514    0.10240 -10.304  < 2e-16 ***
## as.factor(DID)281       -0.82032    0.10204  -8.039 9.05e-16 ***
## as.factor(DID)282       -0.78629    0.10925  -7.197 6.16e-13 ***
## as.factor(DID)284       -0.78193    0.10472  -7.467 8.19e-14 ***
## as.factor(DID)287       -1.04173    0.10478  -9.942  < 2e-16 ***
## as.factor(DID)289       -0.88492    0.10127  -8.738  < 2e-16 ***
## as.factor(DID)290       -1.05728    0.10444 -10.123  < 2e-16 ***
## as.factor(DID)315       -0.74213    0.11562  -6.419 1.37e-10 ***
## as.factor(DID)316       -1.09485    0.10357 -10.571  < 2e-16 ***
## as.factor(DID)318       -1.66751    0.10549 -15.807  < 2e-16 ***
## as.factor(DID)319       -2.28183    0.09897 -23.056  < 2e-16 ***
## as.factor(DID)320       -1.53198    0.10486 -14.609  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 253748  on 198147  degrees of freedom
## Residual deviance: 232377  on 198029  degrees of freedom
## AIC: 232615
## 
## Number of Fisher Scoring iterations: 5
## Warning in ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary =
## TRUE, : maximum three independent variables are allowed
## NULL
Hosmer-Lemeshow
hoslem.test(x = glm_child$y, y = fitted(glm_child))
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  glm_child$y, fitted(glm_child)
## X-squared = 212.13, df = 8, p-value < 2.2e-16

fit the momdel for every districts using apply function

All Districts. C003_01 ~ n_children_in_household PR004_PR009_01 + C002_01

fit_all <- sapply(DID_unique, function(id){
  
  glm_child <- glm(C003_01 ~ n_children_in_household + PR004_PR009_01 + C002_01, 
                   family = "binomial", data = child_ica_dummy %>% filter(DID == id))
  
  # ggPredict(glm_child, se = TRUE, colorAsFactor = TRUE, show.summary = TRUE, point = TRUE)
  
  hoslem <- hoslem.test(x = glm_child$y, y = fitted(glm_child))
  data.frame(id, hoslem["p.value"])
}) %>% t()

fit_all <- fit_all %>% as.data.frame()

fit_all$id <- fit_all$id %>% as.numeric()
fit_all$p.value <- fit_all$p.value %>% as.numeric()

fit_all
Districts whose model is not fitted well
fit_all %>% 
  filter(p.value < 0.05)

memo

fit the model

glm_child <- glm(C003_01~ C001 + C002_01 + PR004_01, family = "binomial", data = child_ica_dummy)

summary(glm_child)
## 
## Call:
## glm(formula = C003_01 ~ C001 + C002_01 + PR004_01, family = "binomial", 
##     data = child_ica_dummy)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.4716  -1.0691   0.6141   0.8525   1.4442  
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.567479   0.012767  -44.45   <2e-16 ***
## C001         0.174066   0.001340  129.90   <2e-16 ***
## C002_01     -0.562952   0.009369  -60.09   <2e-16 ***
## PR004_01     0.788430   0.010620   74.24   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 299726  on 245746  degrees of freedom
## Residual deviance: 271929  on 245743  degrees of freedom
## AIC: 271937
## 
## Number of Fisher Scoring iterations: 4
exponential transformation
exp(glm_child$coefficients)
## (Intercept)        C001     C002_01    PR004_01 
##   0.5669529   1.1901343   0.5695255   2.1999390
confidence interval (intercept and coefficient)
confint(glm_child, level = 0.95)
## Waiting for profiling to be done...
##                  2.5 %     97.5 %
## (Intercept) -0.5925094 -0.5424624
## C001         0.1714426  0.1766953
## C002_01     -0.5813180 -0.5445920
## PR004_01     0.7676346  0.8092649
exponential transformation of confidence interval (odds ratio of intercept and coefficient)
exp(confint(glm_child, level = 0.95))
## Waiting for profiling to be done...
##                 2.5 %    97.5 %
## (Intercept) 0.5529380 0.5813151
## C001        1.1870160 1.1932675
## C002_01     0.5591609 0.5800784
## PR004_01    2.1546637 2.2462561
AIC
extractAIC(glm_child)
## [1]      4.0 271937.4
BIC
extractAIC(glm_child, k = log(nrow(glm_child$data)))
## [1]      4.0 271979.1
effectiveness of explanatory variables
glm_child_null <- glm(C003_01~1, family = "binomial", data = child_ica_dummy)
anova(glm_child_null, glm_child, test = "Chisq")
variables selection
step(glm_child_null, direction = "both", 
     scope = (~ C001 + C002_01 + PR004_01 + PR009_01))
## Start:  AIC=299727.8
## C003_01 ~ 1
## 
##            Df Deviance    AIC
## + C001      1   281139 281143
## + PR009_01  1   293719 293723
## + PR004_01  1   295038 295042
## + C002_01   1   295943 295947
## <none>          299726 299728
## 
## Step:  AIC=281142.5
## C003_01 ~ C001
## 
##            Df Deviance    AIC
## + PR009_01  1   274422 274428
## + PR004_01  1   275560 275566
## + C002_01   1   277803 277809
## <none>          281139 281143
## - C001      1   299726 299728
## 
## Step:  AIC=274428
## C003_01 ~ C001 + PR009_01
## 
##            Df Deviance    AIC
## + C002_01   1   270765 270773
## + PR004_01  1   272824 272832
## <none>          274422 274428
## - PR009_01  1   281139 281143
## - C001      1   293719 293723
## 
## Step:  AIC=270773.3
## C003_01 ~ C001 + PR009_01 + C002_01
## 
##            Df Deviance    AIC
## + PR004_01  1   269072 269082
## <none>          270765 270773
## - C002_01   1   274422 274428
## - PR009_01  1   277803 277809
## - C001      1   289624 289630
## 
## Step:  AIC=269081.6
## C003_01 ~ C001 + PR009_01 + C002_01 + PR004_01
## 
##            Df Deviance    AIC
## <none>          269072 269082
## - PR004_01  1   270765 270773
## - PR009_01  1   271929 271937
## - C002_01   1   272824 272832
## - C001      1   288269 288277
## 
## Call:  glm(formula = C003_01 ~ C001 + PR009_01 + C002_01 + PR004_01, 
##     family = "binomial", data = child_ica_dummy)
## 
## Coefficients:
## (Intercept)         C001     PR009_01      C002_01     PR004_01  
##     -0.7815       0.1760       0.5677      -0.5762       0.4923  
## 
## Degrees of Freedom: 245746 Total (i.e. Null);  245742 Residual
## Null Deviance:       299700 
## Residual Deviance: 269100    AIC: 269100
multicolinearity
vif(glm_child)
##     C001  C002_01 PR004_01 
## 1.010465 1.004067 1.013630